9  BSTR Analysis

This chapter examines reaction engineering of single, isolated, batch stirred-tank reactors (BSTRs). It examines properties, advantages, and disadvantages of BSTRs, common situations where BSTRs are used, the operation of BSTRs, and qualitative and quantitative analysis of BSTRs. It concludes with examples of response and optimization tasks involving BSTRs. Reactor design involving BSTRs is considered separately in Chapter 11.

Simple schematic representation of a BSTR.

9.1 BSTR Characteristics

The two most important characteristics of a BSTR are that (a) the reacting fluid within it is perfectly mixed, and (b) there is no flow of reagents into or out of a BSTR while it is operating, that is, during reaction. As the name batch stirred tank reactor suggests, the most common physical form of a BSTR for chemical processing is a stirred tank. As indicated in the schematic representations in Figure 9.1, they are typically cylindrical in shape, completely containing the reacting fluid (indicated by blue shading) within rigid walls. When processing liquids, there may be some space at the top of the reactor that does not contain reacting fluid. This is referred to as headspace; it is indicated by gray shading in figure. There must be some means of thoroughly and rapidly mixing the reacting fluid. Commonly an agitator is immersed in the reacting fluid for this purpose.

a

b

c
Figure 9.1: Schematic representations of three batch stirred tank reactors, wherein heat is exchanged with an external fluid through the walls of (a) a jacket, (b) a submerged coil, and (c) an external heat exchanger.

If it is necessary to heat or cool the reactor, a heat exchange fluid (indicated by yellow shading in the figure) may flow in and out of the the reactor during operation, but the heat exchange fluid is physically separated from the reacting fluid. Figure 9.1 depicts three common configurations for exchange of heat with an external fluid. One (a in the figure) is referred to as a shell or jacket. It is a compartment that surrounds the reacting fluid and has an external heat exchange fluid flowing through it. Another (b in the figure) consists of a coil of tubing that is submerged in the reacting fluid. The heat exchange fluid flows within the coil. The third (c in the figure) is an external heat exchanger. With the latter configuration it is important that the reacting fluid circulates through the heat exchanger rapidly so that the assumption of perfect mixing of the reacting fluid is satisfied. In all three configurations, it is possible that the exchange fluid temperature is not uniform, but in Reaction Engineering Basics the heat exchange fluid is always assumed to be perfectly mixed.

Stirred tanks are not often used for the processing of gases or when a solid catalyst is being used because stirring gases or small solid particles is difficult. However, configurations different from a stirred tank can be modeled as BSTRs as long as the critical assumptions are met (perfect mixing and no reagent flow in or out of the reactor). Chapter 19 describes a few laboratory BSTRs that are not stirred tanks. Similarly, any reactor in a chemical processing facility that conforms to the assumptions used when deriving the BSTR design equations can be modeled as a BSTR. Heat can be added or removed from a reactor in ways that do not involve an exchange fluid. For example, an electric resistance heating element might be used to add heat. Reaction Engineering Basics only considers heat transfer to or from a perfectly mixed heat exchange fluid.

There are advantages to using a BSTR, as opposed to the other types of ideal reactors, particularly the flow reactors (CSTRs and PFRs). BSTRs offer flexibility both in terms of the reactions being run and in the details of operating the process. In a situation where the demand for products is too small to warrant continuous production, a batch reactor might be used to make one product for a few weeks, and then be used to make a different product for the next few weeks. This kind of flexibility avoids having the reactor sit idle for extended periods of time. Similarly, if the production of a product involves heating to and holding at two or more different temperatures during processing, it may be easier to produce it in a batch reactor than to set up a chain of continuous reactors operating at the necessary temperatures. Similarly, if the addition of one reagent needs to be delayed, it may be possible to do so using a BSTR as described below under BSTR Operation. Again, this may be preferred over two continuous reactors operating in series with that reagent added to the second reactor. In contrast to PFRs, a much broader range of heat transfer area per volume of reacting fluid is possible. Even after the reactor has been constructed, the length of a submerged coil or the area of an external heat exchanger can be altered if additional heat transfer capability is needed.

There are also disadvantages to using BSTRs. As described below, the operation of a BSTR includes periods of time during which no reaction is taking place. Compared to reactors like CSTRs and PFRs that run continuously, this can lead to a lower net rate of production. Additionally, the operation of a BSTR is more labor intensive than operation of a continuous reactor. This adds to the cost of production, meaning a higher product selling price is necessary in order to make a profit. Another issue with batch processes is that of batch-to-batch consistency. As an example, consider the batch production of an expensive perfume. It is critically important that every batch smells the same. Consequently, in addition to being labor intensive, batch processing requires well-trained operators who are careful to precisely follow the specified operational procedure and not deviate from it. Food and beverage production similarly demand high batch-to-batch consistency. This can be even more challenging if, as an example, the process utilizes agricultural reagents that can vary from season to season or year to year.

In light of their advantages and disadvantages, BSTRs are often used to produce value-added products. These are products for which the market demand is smaller, but the price, and therefore the profit per pound, is greater. This is in contrast to commodity products where the demand for large quantities is high, but the price per pound is relatively small. High volume production is usually performed in a continuous reactor (CSTR or PFR).

9.2 BSTR Operation

By its nature, there are times during which reaction is not taking place in a batch reactor. Before reaction can take place, the reactor must first be cleaned (or, in some cases, sterilized) to remove anything “left behind” by the last reaction that was run. Then the reactants must be prepared and charged into the reactor. From this point on, the operation is much like a cooking recipe. A prescribed series of heating and cooling steps is followed, with reaction taking place during every step. A simple example might be that after the reactants are charged, the reactor is heated to a specified temperature, held at that temperature for a specified length of time, and then cooled back to room temperature. In almost all situations, the reaction has ceased by the end of the last step, either because one of the reactants has been fully consumed or because the temperature has been lowered to the point where the rate is effectively zero. The last step in the operating protocol then is to remove the product from the reactor. Often it is transferred to a holding tank from which it can be sent elsewhere for further processing, purification, or longer term storage.

Of course, the operational protocol, i. e. the “recipe,” could be more complicated, requiring any number of steps (also called stages). In some cases an intermediate stage in the operational protocol might call for the addition of a reagent. As long as that reagent is added instantaneously, the reactor can be analyzed as a BSTR. However, if the reagent is added over time as reactions continue to occur, the reactor is no longer a BSTR and that phase of processing must be analyzed as a semi-batch reactor (Chapter 10).

It is very important to recognize that the productivity of the reactor, that is the amount of product it produces per hour, must be calculated taking the so-called turnaround time into account. That is, one must include the time needed to clean, fill and drain the reactor, not just the length of time during which reaction takes place. Doing so permits the definition of the net rate of generation of reagent \(i\) given in Equation 9.1.

\[ r_{i,net} = \frac{n_{i,f} - n_{i,0}}{t_{rxn} + t_{turn}} \tag{9.1}\]

Even if the turnaround time was zero, it would still be desirable to define a net rate as above. The net rate is useful for characterizing the overall process. As the reaction is progressing, the rate of reaction is continually changing because the temperature and composition are changing. The rate of reaction at any one time during the processing is referred to as the instantaneous rate at that time. The net rate defined in Equation 9.1 is a kind of average of the instantaneous rate over the processing time. Similarly, it is sometimes useful to differentiate between overall and instantaneous values of other quantities such as yields, conversions, selectivities, etc.

9.3 Qualitative Analysis of Reaction in a BSTR

A simple process for the qualitative analysis of a reacting system was presented in Chapter 7.1. In a BSTR, there are no spatial variations in temperature or composition because the reactor is perfectly mixed. Therefore, qualitative analysis of a BSTR is used to determine how temperature, reaction rate, composition and related quantities (conversion, selectivity) vary with reaction time.

For purposes of illustration, suppose the typical, irreversible, exothermic reaction, \(A \rightarrow Z\) takes place in an isothermal BSTR. In this situation, only the concentrations need to be considered because in an isothermal reactor the temperature does not change. At the start of the process, the concentration of A will be high and that of Z will be small or zero. During a very small interval of time at the start of the process the concentration of A will decrease and the concentration of Z will increase due to the reaction starting to take place. Hence, in a plot of the concentration of reactant A versus reaction time, the \(y\)-intercept will be high and the slope at that point will be negative. For the product, Z, the \(y\)-intercept will be low or zero and the initial slope will be positive.

The changes in composition will affect the reaction rate. For a typical, irreversible reaction, the decrease in the reactant concentration will tend to decrease the rate while the increase in the product concentration will not affect it. Consequently, at the end of that small interval of reaction time, the rate will be smaller.

In the next small interval of time, the concentration of the reactant, A, will decrease by less than it did in the first interval because the rate is lower. For the same reason, the concentration of the product, Z, will increase by less than it did during the first interval of time. In terms of graphs showing the concentrations versus time, the concentration of the reactant, A, will start high with a negative slope and a curvature that is concave upward. That is, the slope will become less steep as the reaction time increases. In like manner, a graph showing the concentration of Z versus reaction time will start at or near zero with a positive slope and a curvature that is concave downward (i. e. the slope again will get less steep as the reaction time increases).

If that initial behavior continued throughout the process, the concentration of the reactant would decrease with a smaller and smaller slope until finally reaching a horizontal plateau at a value of zero (because the reaction is irreversible and therefore goes to completion). The concentration of the product would increase with a smaller and smaller slope until finally reaching a plateau at a value equal to the initial concentration of A (because in the reaction being considered, the stoichiometry is 1 Z per A). The results of this qualitative analysis are shown in Figure 9.2. Qualitative analysis is not limited to reagent concentrations. Graphs showing conversion, selectivity, rate, temperature, etc. vs. time can also be generated.

Figure 9.2: Qualitative reactant (left) and product (right) concentration profiles during isothermal reaction in a BSTR.

The qualitative analysis here was straightforward because the reactor was isothermal. In Example 9.6.1 the qualitative analysis of an adiabatic BSTR is presented. That analysis is more complicated because both concentrations and temperature are changing. In fact, the changes in composition and temperature, taken separately, would affect the rate in opposite ways. As a consequence there is greater uncertainty in the qualitative analysis results.

9.4 BSTR Design Equations

The reactor design equations for BSTRs are derived in Appendix E.2 and discussed in Chapter 6.4. The BSTR mole and energy balances shown below were presented in Equations 6.8 and 6.9. Chapter 6.3 also presented the energy balances shown below for a heat exchange fluid that transfers sensible heat, Equation 6.2, and one that transfers latent heat, Equation 6.6. That chapter also describes how to select the specific reactor design equations needed to model BSTRs under different circumstances. That discussion will not be duplicated here.

\[ \frac{dn_i}{dt} = V \sum_j \nu_{i,j}r_j \]

\[ \left(\sum_i n_i \hat C_{p,i} \right) \frac{dT}{dt} - V\frac{dP}{dt} - P \frac{dV}{dt} = \dot Q - \dot W - V \sum_j \left(r_j \Delta H_j \right) \]

\[ \rho_{ex} V_{ex} \tilde C_{p,ex}\frac{dT_{ex,out}}{dt} = -\dot Q - \dot m_{ex} \int_{T_{ex,in}}^{T_{ex,out}} \tilde C_{p,ex}dT \]

\[ \frac{\rho_{ex} V_{ex} \Delta H_{\text{latent},ex}^0}{M_{ex}} \frac{d \gamma}{dt} = - \dot Q - \gamma \dot m_{ex} \frac{\Delta H_{\text{latent},ex}^0}{M_{ex}} \]

Note

It is important to notice that the \(V\frac{dP}{dt}\) and \(P\frac{dV}{dt}\) terms in the reacting fluid energy balance will have units of pressure volume time-1. Most of the other terms will have units of energy time-1. This means that a unit conversion is needed. One easy way to do this is to multiply the \(V\frac{dP}{dt}\) and \(P\frac{dV}{dt}\) terms by the ideal gas constant in units of energy mol-1 temperature-1 and divide by the ideal gas constant in units of pressure volume mol-1 temperature-1. The moles and temperatures cancel out giving the necessary conversion factor. This works even for a liquid phase system because the ideal gas law is not being used as the equation of state (i. e. to convert from pressure volume to mol temperature). In Reaction Engineering Basics the rate of doing work, \(\dot{W}\) is almost always negligible. However, if it is non-zero, it will likely have units of power (e. g. hp) or pressure volume time-1, and again a unit conversion is needed.

The design equations can often be simplified. As noted in Chapter 6.4, if the reactor walls are rigid (or if an ideal, incompressible liquid is being processed) then the volume, \(V\), will be constant and its time derivative, \(\frac{dV}{dt}\), will equal zero. If the work associated with agitation is negligible (it usually is) and there are no other shafts or moving boundaries, then the rate at which the reacting fluid performs work, \(\dot W\), on its surroundings is also equal to zero. When an ideal, incompressible liquid is being processed, the pressure, \(P\), is constant and its time derivative, \(\frac{dP}{dt}\), is equal to zero. (The pressure gradient within the reacting fluid due to the hydrostatic head is neglibible.) Finally, for liquids, Equation 6.10 showed how the sensible heat term in Equation 6.9 can be expressed in terms of the mass-specific or volume-specific heat capacity of the solution as a whole instead of the individual molar heat capacities.

\[ \left(\sum_i n_i \hat C_{p,i} \right) \frac{dT}{dt}\ \Leftrightarrow\ \rho V \tilde C_p \frac{dT}{dt}\ \Leftrightarrow\ V \breve C_p \frac{dT}{dt} \]

When rate expressions are substituted into the BSTR mole and energy balances, they typically introduce concentrations or partial pressures of reagents. When the design equations are solved, those concentrations and partial pressures must be expressed in terms of the molar amounts of the reagents (see Appendix F.5.2). The concentration is simply the molar amount divided by the volume, Equation 2.7.

\[ C_i = \frac{n_i}{V} \]

For ideal gases, the partial pressure is related to the molar amount through the ideal gas law, Equation 2.12.

\[ P_i = \frac{n_iRT}{V} \]

When heat transfer involves an exchange fluid as shown in Figure 9.1, the rate of heat exchange in the energy balances can be calculated using Equation 6.18. If the reactor operates adiabatically (i. e. it is perfectly insulated and does not have a shell, coil or external heat exchanger), \(\dot{Q}\) is equal to zero.

\[ \dot Q = UA\left( T_{ex,out} - T \right) \]

9.5 General Approach to Modeling Isolated BSTRs

Chapter 7 described how to identify a reaction engineering assignment that involves the modeling of an isolated ideal reactor, and Chapter 8 suggested a six step workflow for reactor analysis.

  1. Summarize the information provided in the assignment.
  2. Generate the model equations.
  3. Generate equations for calculating the quantities of interest.
  4. Formulate the calculations.
  5. Implement the calculations numerically and execute them.
  6. Present and discuss the results.

This section describes the application of that workflow to the modeling of an isolated BSTR. The process may seem a bit abstract when first read, but it should become clearer after seeing it used in the examples provided at the end of this chapter.

Often the operational protocol for a BSTR has multiple heating and cooling stages. When this is true, each stage of the protocol may need to be analyzed separately from the other stages because terms in the reactor design equations may change. The reactor design equations for BSTRs are ordinary differential equations with elapsed time, \(t\), as the independent variable. The stages occur sequentially, and consequently the values of \(t\) and the dependent variables at the end of one stage in the operating protocol become the initial values for the next stage. The one exception to this is if a reagent is instantaneously added at the end of a stage. In that situation, the initial values for the next stage are the result of mixing the added reagent with the reagents in the BSTR at the end of the previous stage.

Summarizing the information provided in the assignment. An assignment summary is useful in that it extracts the critical information from the assignment narrative and presents it as a concise summary. Assignment summaries are described in Chapter 8.1. When analyzing an isolated BSTR the assignment summary will typically include lists of the reactions taking place, the rate expressions for those reactions, the quantities of interest, known constants given in the assignment narrative, and any specified parameters. The summary will also include details of the operational protocol for the BSTR, a schematic diagram of it, and the basis for the calculations, if there is one.

Generating the model equations. When analyzing an isolated BSTR , it will always be necessary to generate the reactor design equations that are needed to model that specific BSTR. Selecting which equations are needed to model a given BSTR and simplifying the design equations is described in Chapter 6.4 and in the preceding section of this chapter.

The design equations for an isolated BSTR will be either a set of IVODEs or a set of DAEs, with \(t\) as the independent variable. The dependent variables will always include the molar amounts of the reagents present in the system, \(n_i\). Depending which other design equation are needed to model the BSTR, the dependent variables may also include the temperature of the reacting fluid, \(T\), its pressure, \(P\), and either the instantaneous outlet temperature of the exchange fluid, \(T_{ex,out}\), or the instantaneous fraction of the exchange fluid undergoing phase change, \(\gamma\).

If the reacting fluid is an ideal gas and the reactor design equations include an energy balance on the reacting fluid, the temperature and the pressure will appear as a dependent variables. Consequently, the number of dependent variables will be one greater than the number of IVODEs. In this situation either a dependent variable must be eliminated from the design equations or another IVODE must be added to them (Appendix F.5.2). A differential form of the ideal gas law, Equation 6.11, is used to accomplish that. For a BSTR, \(V\) is constant so its time derivative is zero, and does not appear in the differential form of the gas law. The differential ideal gas law can either be added to the reactor design equations, making the number of equations equal to the number of dependent variables, or it can be solved for \(\frac{dP}{dt}\) and substituted into the energy balance, thereby eliminating \(P\) as a dependent variable.

\[ RT\left( \sum_i \frac{dn_i}{dt} \right) + R\left( \sum_i n_i \right)\frac{dT}{dt} - V \frac{dP}{dt} = 0 \]

The reactor design equations will be solved numerically, so if they are not in the form of derivative expressions, they should be converted to that form (see Appendix F.4.2).

Next it is useful to create a table that identifies the initial values and known final values of the ODE variables. For the first stage of the operational protocol, the initial value of \(t\) can usually be set equal to 0.0. The other initial values and the known final values should be indicated using the same variable symbol as that used in the assignment summary. Quite commonly, a subscripted zero is added to indicate an initial value and a subscripted “f” is used to denote a final value.

The table can then be used to write supporting equations to calculate values that are not known. If only one final value is known or can be calculated directly, the design equations are IVODEs. If two final values are known or can be calculated directly, the design equations are DAEs and an implicit equation can be written relating one of the final values to either an unknown initial value or an unknown constant appearing in the ODEs (see Appendix F.5.3).

Finally, any ancillary equations that are needed to solve the BSTR design equations can be written. These will include rate expressions, equations for calculating the rate coefficients and concentrations appearing in them, and equations for calculating any other additional unknowns that appear in the design equations using the independent and dependent variables and known constants.

Calculating the quantities of interest. The BSTR design equations will be solved numerically, yielding corresponding sets of values of the independent and dependent variables spanning the range from their initial to their final values. These may not be the quantities of interest in the assignment, but they can be used to calculate the quantities of interest. Accordingly, any equations for calculating the quantities of interest should be generated.

Formulating the calculations. Once all of the necessary equations for completing the assignment have been generated, the calculations can be formulated. Basically, as noted in Chapter 8.4, this entails writing specifications for computer functions that will solve the equations numerically. For an isolated BSTR the design equations will either be IVODEs or DAEs. Writing specifications for the reactor function, derivatives function, and residuals function that are needed for solving the design equations is described in Appendices F.5.2 and F.5.3. Finally, specifications should be written for a “deliverables” function that calls those functions to solve the design equations and then uses the results to calculate the quantities of interest.

Completing the assignment workflow. The computer function specifications included in the formulation of the calculations provide essentially all of the information that is needed for implementing and performing them. The choice of the specific software package and programming language for doing so is left to the reader. Examples in Reaction Engineering Basics will simply note that the calculations were performed as described in the formulation and proceed to presenting and discussing the results. As noted in Chapter 8.6 a qualitative analysis of the BSTR can provide insights when discussing the results.

9.6 Examples

The examples in this section illustrate the analysis of isolated BSTRs. Examples 9.6.1 and 9.6.2 are response assignments where the assignment specifies a sufficient number of quantities so that the reactor design equations can be solved to find the quantities of interest. In Example 9.6.1 the BSTR design equations are IVODEs while in Example 9.6.2 they are DAEs. In addition to illustrating the computational procedure for these two types of equations, Example 9.6.1 involves a single reaction taking place in an adiabatic BSTR while Example 9.6.2 features two reactions taking place in a BSTR that is cooled using an exchange fluid.

Examples 9.6.3 and 9.6.4 are optimization assignments that identify a reactor variable that must be maximized or minimized by specifying a reactor input. The BSTR operational protocol in Example 9.6.3, like those in Examples 9.6.1 and 9.6.2, consists of a single stage, but unlike the other examples, the reacting fluid is a gas and not a liquid. Example 9.6.4 illustrates a situation where there are multiple stages in the operating protocol. It also considers turnaround time and the net rate of reaction.

9.6.1 Concentration, Temperature and Rate Profiles during Adiabatic Operation of a BSTR

The liquid-phase reaction between A and B, equation (1), is irreversible. At the conditions of interest, the heat of reaction is constant and equal to -101.2 kJ mol-1. The rate expression is given in equation (2), where the rate coefficient exhibits Arrhenius temperature dependence with a pre-exponential factor of 5.11 x 104 L mol-1 s-1 and an activation energy of 74.8 kJ mol-1.

Two solutions, one containing only reagent A at 180 °C and one containing only reagent B at 180 °C, are charged to an adiabatic BSTR producing 1900 L of solution initially containing 2.9 mol A L-1 and 3.2 mol B L-1 at 180 °C. The solution is ideal and has a constant heat capacity of 1.23 cal g-1 K-1 and a constant density of 1.02 g cm-3. Plot the concentrations of A, B, Y, and Z, the reacting fluid temperature, and the instantaneous reaction rate for the first 2 h of reaction, and comment upon the shapes of the graphs.

\[ A + B \rightarrow Y + Z \tag{1} \]

\[ r = kC_AC_B \tag{2} \]


This assignment describes a single BSTR and no other equipment. Information about the reaction, rate, and reactor operation are provided. The quantities of interest are related to the final reagent amounts and temperature. This indicates I will need to solve the BSTR design equations for the reactor.

I’ll start by summarizing the information provided in the assignment narrative. I will use a subscripted “0” to denote initial values and a subscripted “f” to denote final values. I’ll also add the ideal gas constant to the given constants. Numerical solution of the design equations will yield sets of corresponding values of the independent and dependent variables that span the range from their initial values to their final values. I’ll use underlined variable symbols, e. g. \(\underline{t}\), to represent sets of values.

9.6.1.1 Assignment Summary

Reaction:

\[ A + B \rightarrow Y + Z \tag{1} \]

Rate Expression:

\[ r = kC_AC_B \tag{2} \]

Reactor System: Adiabatic BSTR

Reactor Schematic:

Figure 9.3: A closed (batch), adiabatic stirred tank reactor.

Quantities of Interest: \(\underline{C}_A\), \(\underline{C}_B\), \(\underline{C}_Y\), \(\underline{C}_Z\), \(\underline{T}\), and \(\underline{r}\) vs. \(\underline{t}\) as graphs.

Given and Known Constants: \(\Delta H_1\) = -101.2 kJ mol-1, \(k_{0}\) = 5.11 x 104 L mol-1 s-1, \(E\) = 74.8 kJ mol-1, \(V\) = 1900 L, \(T_0\) = (180 + 273.15) K, \(C_{A,0}\) = 2.9 mol L-1, \(C_{B,0}\) = 3.2 mol L-1, \(\tilde{C}_p\) = 1.23 cal g-1 K-1, \(\rho\) = 1.02 g cm-3, \(t_f\) = 2 h, and \(R\) = 8.314 J mol-1 K-1.

9.6.1.2 Generation of the Model Equations

Reactor analysis always requires a set of reactor design equations. The selection of the specific design equations for a given analysis was described in Chapter 6.4. Mole balances are always needed. I’ll write a BSTR mole balance, Equation 6.8, on each of the reagents. In this assignment, only one reaction is taking place, so the summation reduces to a single term.

\[ \frac{dn_i}{dt} = V \sum_j \nu_{i,j}r_j \qquad \Rightarrow \qquad \frac{dn_i}{dt} = V \nu_i r \]

This BSTR is not isothermal, so the mole balances cannot be solved independently of an energy balance on the reacting fluid, Equation 6.9. In this assignment the reacting fluid is a liquid which I will assume to be an incompressible ideal solution. Consequently, the reacting fluid volume and the total pressure are constant, and their time derivatives are equal to zero. The reactor is adiabatic, so the rate of heat input, \(\dot{Q}\), is equal to zero. There are no shafts or moving boundaries other than the agitator, and the rate at which it performs work is negligible, so \(\dot{W}\) is also equal to zero. With only one reaction taking place, the summation over the reactions reduces to a single term. This assignment provides the mass-specific heat capacity of the entire solution, so the sensible heat term can be re-written using that heat capacity, Equation 6.10. Finally, the energy balance can be put in the form of a derivative expression by dividing both sides by \(\rho V \tilde{C}_p\). Since only one reaction is taking place, I don’t need a subscript on the rate.

\[ \cancelto{\rho V \tilde{C}_p}{\left(\sum_i n_i \hat C_{p,i} \right)} \frac{dT}{dt} - \cancelto{0}{V\frac{dP}{dt}} - \cancelto{0}{P\frac{dV}{dt}} = \cancelto{0}{\dot Q} - \cancelto{0}{\dot W} - V \cancelto{r \Delta H}{\sum_j \left(r_j \Delta H_j \right)} \]

\[ \rho V \tilde{C}_p \frac{dT}{dt} = - V r \Delta H_1 \]

\[ \frac{dT}{dt} = -\frac{V r \Delta H_1}{\rho V \tilde{C}_p} \]

In this assignment the BSTR operates adiabatically, so there isn’t an exchange fluid, and an exchange fluid energy balance cannot be written. Momentum balances are not used with stirred tank reactors. This gives five reactor design equations that contain five dependent variables, and consequently it is not necessary to add an IVODE or eliminate a dependent variable. So for this reactor, the design equation consist of mole balances on A, B, Y, and Z and an energy balance on the reacting fluid.

BSTR Design Equations

Mole balance design equations for A, B, Y, and Z are presented in equations (3) through (6), and the energy balance on the reacting fluid is given by equation (7).

\[ \frac{dn_A}{dt} = -rV \tag{3} \]

\[ \frac{dn_B}{dt} = -rV \tag{4} \]

\[ \frac{dn_Y}{dt} = rV \tag{5} \]

\[ \frac{dn_Z}{dt} = rV \tag{6} \]

\[ \frac{dT}{dt} = -\frac{r \Delta H_1}{\rho \tilde{C}_p} \tag{7} \]

Initial values and a stopping criterion are needed when solving IVODEs or DAEs. The instant at which the two solutions are charged into the reactor can be defined as \(t=0\). The initial values are then the molar amounts of the reagents and the temperature at that instant. The two solutions that were charged to the BSTR did not contain reagent Y or reagent Z, so their initial amounts are zero. The final value of the time in this assignment is 2 h, which I represented as \(t_f\) in the list of known constants. Consequently the stopping criterion is that \(t=t_f\).

No other final values are known or can be calculated directly, so the design equations are IVODEs. The initial temperature is known, but the initial molar amounts of A and B are not. However they can be calculated directly from the initial concentrations and the fluid volume, which are known.

Intial Values and Stopping Criterion

Table 9.1: Initial values and stopping criterion for solving the design equations, equations (3) through (7).
Variable Initial Value Stopping Criterion
\(t\) \(0\) \(t_f\)
\(n_A\) \(n_{A,0}\)
\(n_B\) \(n_{B,0}\)
\(n_Y\) \(0\)
\(n_Z\) \(0\)
\(T\) \(T_0\)

Supporting Equations for Calculating the Initial and Final Values

The initial molar amounts of A and B in Table 9.1 can be calculated using equations (8) and (9).

\[ n_{A,0} = C_{A,0}V \tag{8} \]

\[ n_{B,0} = C_{B,0}V \tag{9} \]

To solve the IVODEs numerically I’ll need to write a derivatives function. It will receive \(t\), \(n_A\), \(n_B\), \(n_Y\), \(n_Z\), and \(T\) for the current integration step as arguments, and must return the derivatives calculated using equations (3) through (7). Looking at those equations, all of the other quantities except \(r\) are known, so before I can evaluate the derivatives, I need to calculate \(r\).

The assignment provides the rate expression, equation (2), which can be used to calculate \(r\), but the rate expression contains the rate coefficient, \(k\), and the concentrations of A and B. They, too, will need to be calculated before \(r\) can be calculated. The rate coefficient can be calculated using the Arrhenius expression, Equation 4.8. The defining equation for concentration in a closed system, Equation 2.7, can be used to calculate the concentrations of A and B.

Ancillary Equations for Evaluating the Derivatives

The rate, \(r\), appearing in the reactor design equations can be calculated using the given rate expression, equation (2). The rate coefficient, \(k\), appearing in equation (2) can be calculated using the Arrhenius expression, equation (10). The concentrations can be calculated using equations (11) and (12).

\[ k = k_{0}\exp{\left( \frac{-E}{RT} \right)} \tag{10} \]

\[ C_A = \frac{n_A}{V} \tag{11} \]

\[ C_B = \frac{n_B}{V} \tag{12} \]

At this point the IVODEs can be solved numerically to obtain corresponding sets of values of \(t\), \(n_A\), \(n_B\), \(n_Y\), \(n_Z\), and \(T\) spanning the range from their initial to their final values. \(\underline{T}\) and \(\underline{t}\) are two of the quantities of interest needed to make one of the requested graphs, but the others need to be calculated using the results from solving the IVODEs. Again, the defining equation for concentration in a closed system, Equation 2.7, can be used to calculate the concentrations. The instantaneous rate can then be calculated using the rate expression.

9.6.1.3 Calculating the Quantities of Interest

With the information given above, equations (3) through (7) can be solved for corresponding sets of values of \(t\), \(n_A\), \(n_B\), \(n_Y\), \(n_Z\), and \(T\) that span the range from their initial values to their final values. The concentrations of A, B, Y, and Z, at each time can then be calculated using equations (11), (12), (13), and (14). Then, knowing \(C_A\), \(C_B\), and \(T\), at each time, the instantaeous rate in that range can be calculated using the Arrhenius expression, equation (10), and the rate expression, equation (2).

\[ C_Y = \frac{n_Y}{V} \tag{13} \]

\[ C_Z = \frac{n_B}{V} \tag{14} \]

9.6.1.4 Formulation of the Calculations

To complete the assignment, the BSTR design equations are solved first, and then the quantities of interest are calculated and the requested graphs are generated. The BSTR design equations are IVODEs, and they can be solved numerically as described in Appendix F.5.2 by writing a derivatives function and a BSTR function. A deliverables function can be added to calculate the quantities of interest and generate the graphs.

  • The derivatives function must receive the independent and dependent variables at the start of an integration step as arguments, and it must return the derivatives of the dependent variables. To do so it must calculate the rate before it can evaluate the derivatives.
  • The BSTR function does not need any arguments to calculate and return corresponding sets of values of the independent and dependent variables spanning the interval from \(t_0\) to \(t_f\). To do so it must define the initial values and stopping criterion and then pass them to an IVODE solver along with the derivatives function.
  • The deliverables function does not need any arguments, and it doesn’t return anything. It calls the BSTR function, uses the results it returns to calculate the quantities of interest, and generates the requested graphs.

Specifications can be written for each of these functions in preparation for writing them. The function specifications given here assume that the known constants are available wherever they are needed.

Derivatives Function Specifications

Arguments: \(t\), \(n_A\), \(n_B\), \(n_Y\), \(n_Z\), and \(T\)

Returns: \(\frac{dn_A}{dt}\), \(\frac{dn_B}{dt}\), \(\frac{dn_Y}{dt}\), \(\frac{dn_Z}{dt}\), and \(\frac{dT}{dt}\)

Algorithm:

\[ k = k_{0}\exp{\left( \frac{-E}{RT} \right)} \tag{10} \]

\[ C_A = \frac{n_A}{V} \tag{11} \]

\[ C_B = \frac{n_B}{V} \tag{12} \]

\[ r = kC_AC_B \tag{2} \]

\[ \frac{dn_A}{dt} = -rV \tag{3} \]

\[ \frac{dn_B}{dt} = -rV \tag{4} \]

\[ \frac{dn_Y}{dt} = rV \tag{5} \]

\[ \frac{dn_Z}{dt} = rV \tag{6} \]

\[ \frac{dT}{dt} = -\frac{r \Delta H_1}{\rho \tilde{C}_p} \tag{7} \]

BSTR Function Specifications

Arguments: none

Returns: \(\underline{t}\), \(\underline{n}_A\), \(\underline{n}_B\), \(\underline{n}_Y\), \(\underline{n}_Z\), and \(\underline{T}\)

Algorithm:

\[ t_{initial} =0 \tag{15} \]

\[ n_{A,initial} = n_{A,0} = C_{A,0}V \tag{8} \]

\[ n_{B,initial} = n_{B,0} = C_{B,0}V \tag{9} \]

\[ n_{Y,initial} = 0 \tag{16} \]

\[ n_{Z,initial} = 0 \tag{17} \]

\[ T_{initial} = T_0 \tag{18} \]

\[ t_{final} = t_f \tag{19} \]

 

\[ \begin{matrix} t_{initial}, n_{A,initial}, n_{B,initial}, n_{Y,initial}, \\ n_{Z,initial}, T_{initial}, t_{final}, \text{derivatives function}\\ \Downarrow \\ \text{IVODE solver} \\ \Downarrow \\ \underline{t}, \underline{n}_A, \underline{n}_B, \underline{n}_Y, \underline{n}_Z, \underline{T} \end{matrix} \tag{20} \]

Deliverables Function Specifications

Arguments: none

Returns: nothing

Algorithm:

\[ \begin{matrix} \text{BSTR function} \\ \Downarrow \\ \underline{t}, \underline{n}_A, \underline{n}_B, \underline{n}_Y, \underline{n}_Z, \underline{T} \end{matrix} \tag{21} \]

 

\[ \underline{C}_A = \frac{\underline{n}_A}{V} \tag{11} \]

\[ \underline{C}_B = \frac{\underline{n}_B}{V} \tag{12} \]

\[ \underline{C}_Y = \frac{\underline{n}_Y}{V} \tag{13} \]

\[ \underline{C}_Z = \frac{\underline{n}_Z}{V} \tag{14} \]

\[ \underline{k} = k_{0}\exp{\left( \frac{-E}{R\underline{T}} \right)} \tag{10} \]

\[ \underline{r} = \underline{k} \underline{C}_A \underline{C}_B \tag{2} \]

 

\[ \begin{matrix} \underline{t}, \underline{C}_A, \underline{C}_B, \\ \underline{C}_Y, \underline{C}_Z, \underline{T} \\ \Downarrow \\ \text{plotting function} \\ \Downarrow \\ \text{requested graphs} \end{matrix} \tag{22} \]

9.6.1.5 Implementation and Execution of the Calculations

The computer functions above were written as specified using an appropriate mathematics software package. The calculations were then performed by calling the deliverables function.

9.6.1.6 Results and Discussion

The variation of concentrations of reagents, A, B, Y, and Z, during the first two hours of operation of the BSTR are shown in Figure 9.4. The variation of the reacting fluid temperature during that period is shown in Figure 9.5, and the variation of the instantaneous rate is shown in Figure 9.6.

Figure 9.4: Concentrations versus reaction time.
Figure 9.5: Temperature versus reaction time.
Figure 9.6: Instantaneous rate versus reaction time.

The concentrations of the reactants, A and B, decrease monotonically during the first two hours of operation while the concentrations of the products, Y and Z, increase. The temperature increases steadily, too. The variation of the reaction rate is different. Initially it increases, but then it passes through a maximum after which it steadily decreases. A qualitative analysis shows that the variations in the concentrations, temperature and rate seen in Figures 9.4, 9.5, and 9.6 are expected. At \(t=0\), the concentrations of the reagents and the temperature are all at their initial values.

At that point the rate is positive, so after a short interval of time the concentrations of A and B will have decreased slightly and the concentrations of Y and Z will have increased slightly due to the occurrence of the reaction. Because the reaction is exothermic and the reactor operates adiabatically, the temperature will have increased due to the heat released by the reaction.

The rate could either increase or decrease during this short interval. The concentrations of the reactants, A and B, are decreasing, and this alone would cause the rate to decrease. However, at the same time the temperature is increasing, and that alone would cause the rate to increase. Generally, at the start of the reaction the exponential dependence of the rate upon the temperature will be stronger than its dependence upon the reactant concentrations, and the rate will increase. This is seen to be the case in Figure 9.6. (Had the heat of reaction been very small, leading to a small increase in the temperature, it is possible that the rate would have decreased, especially if the activation was also small).

To summarize, during that initial short interval of time, the concentrations of Y and Z, the temperature, and the rate all have positive slopes while the concentrations of A and B have negative slopes. During the next short interval of time, the rate is larger than during the first interval of time. As a consequence, the concentrations and temperature will change more during the second interval. That means that the curvature of the temperature, Y concentration, and Z concentration profiles is initially concave upward and the curvature of the A and B concentration profiles is initially concave downward. Again during this second short interval, opposing effects on the rate are associated with the decreasing reactant concentrations and the increasing temperature.

It is hard to see in Figure 9.4 and Figure 9.5, but because initially the rate is increasing, the concentrations of A and B are decreasing faster and faster and the concentrations of Y and Z and the temperature are increasing faster and faster. These trends cannot continue indefinitely, if they did, \(C_Y\), \(C_Z\), and \(T\) would approach positive infinity and \(C_A\) and \(C_B\) would approach negative infinity.

In fact, the curvature does change because at some point the effect of decreasing reactant concentrations upon the rate becomes stronger than the effect of increasing temperature. At the point where the two effects become equal, the rate reaches a maximum. Beyond the maximum, the effect of decreasing reactant concentration dominates and the rate decreases. This is seen in Figure 9.6 where the instantaneous rate passes through a maximum.

At the point where the rate reaches its maximum value, the concentration and temperature profiles exhibit an inflection point. (Again, in this example this is very, very subtle and cannot be discerned in the figures.) Beyond that point, the \(C_Y\), \(C_Z\), and \(T\) profiles are concave downward and the \(C_A\) and \(C_B\) are concave upward. That is, the changes over time become smaller and smaller. This behavior can continue indefinitely. Eventually the concentrations and temperature will stop changing and become constant. In this example, A is the limiting reactant, and the reaction is irreversible, so the concentrations and temperature will become constant when \(C_A\) becomes equal to zero.

9.6.2 Parallel Reaction Selectivity in a Cooled BSTR with an Unknown Initial Temperature

A jacketed, 10 L BSTR is going to be used to process an aqueous solution of A and B. The reactor jacket is perfectly mixed with a volume of 1400 cm3, a heat transfer coefficient of 138 cal ft−2 min−1 K−1 and a heat transfer area of 1200 cm2. Cooling water at 40 °C flows into the jacket at a rate of 100 g min-1. Both the reacting fluid and the cooling water may be taken to have a density of 1 g cm−3 and a heat capacity of 1 cal g−1 K−1.

Solutions containing A, B, X, Y, and Z are ideal and there is no heat of mixing. The reactor will be charged by mixing two separate solutions, one containing only reagent A and the other containing only reagent B, giving a 10 L charge with initial concentrations of 5 mol A L-1 and 7 mol B L-1. At the time the reactor is charged, the cooling water temperature is 40 °C. Reactions (1) and (2) will occur in the reactor. The heat of reaction (1) is −16.7 kcal mol−1, and that for reaction (2) is −14.3 kcal mol−1. The rate expression for reaction (1) is given in equation (3) where the pre-exponential factor is 9.74 x 109 L mol−1 min−1 and the activation energy is 20.1 kcal mol−1. The rate expression for reaction (2) is given in equation (4) where the pre-exponential factor is 2.38 x 1013 min−1 and the activation energy is 25.3 kcal mol−1.

What initial temperature is needed in order to convert 45% of the A in 30 min, and with that initial temperature what will the final temperature, the outlet exchange fluid temperature, and the selectivity for X over Z equal? Discuss ways to increase the selectivity.

\[ A + B \rightarrow X + Y \tag{1} \]

\[ A \rightarrow Z \tag{2} \]

\[ r_1 = k_1C_AC_B \tag{3} \]

\[ r_2 = k_2C_A \tag{4} \]


This is an isolated reactor modeling assignment: the system consists of only a reactor, it’s operation is described, the rate expressions are known and it asks about reactor inputs and outputs. I’ll begin by summarizing the assignment. I’ll use appropriate variable symbols to represent the given constants. For initial values I’ll add a subscripted “0” and for final values a subscripted “f”. Cooling water flows into the jacket continuously, so I’ll denote the temperature of the cooling water entering the jacket as \(T_{ex,in}\) and the temperature leaving the jacket as \(T_{ex,out}\). Because the jacket is perfectly mixed, \(T_{ex,out,0}\) is the initial temperature of the cooling water in the jacket. I will add the ideal gas constant to the given constants.

9.6.2.1 Assignment Summary

Reactions:

\[ A + B \rightarrow X + Y \tag{1} \]

\[ A \rightarrow Z \tag{2} \]

Rate Expressions:

\[ r_1 = k_1C_AC_B \tag{3} \]

\[ r_2 = k_2C_A \tag{4} \]

Reactor System: Cooled BSTR.

Reactor Schematic:

Figure 9.7: A closed (batch), non-adiabatic stirred tank reactor.

Quantities of Interest: \(T_0\), \(T_f\), \(T_{ex,f}\), and \(S_{X/Z,f}\)

Given and Known Constants: \(V\) = 10 L, \(V_{ex}\) = 1.4 L, \(U\) = 138 cal ft-2 min-1 K-1, \(A\) = 1200 cm2, \(T_{ex,in}\) = 40 °C, \(\dot{m}_{ex}\) = 100 g min-1, \(\rho\) = 1.0 g cm-3, \(\rho_{ex}\) = 1.0 g cm-3, \(\tilde{C}_p\) = 1.0 cal g-1 K-1, \(\tilde{C}_{p,ex}\) = 1.0 cal g-1 K-1, \(C_{A,0}\) = 5.0 M, \(C_{B,0}\) = 7.0 M, \(T_{ex,out,0}\) = 40 °C, \(\Delta H_1\) = -16.7 kcal mol-1, \(\Delta H_2\) = -14.3 kcal mol-1, \(k_{0,1}\) = 9.74 x 109 L mol-1 min-1, \(E_1\) = 20.1 kcal mol-1, \(k_{0,2}\) = 2.38 x 1013 min-1, \(E_2\) = 25.3 kcal mol-1, \(f_{A,f}\) = 0.45, \(t_f\) = 30 min, and \(R\) = 1.987 cal mol-1 K-1.

9.6.2.2 Generation of the Model Equations

The first thing I need to do is generate the reactor design equations for this reactor. Mole balances are always needed, and for a BSTR the mole balance is Equation 6.8. In this system, there are two reactions taking place, so the summation expands to two terms.

\[ \frac{dn_i}{dt} = V \sum_j \nu_{i,j}r_j = \nu_{i,1}r_1V + \nu_{i,2}r_2V \]

The BSTR is not isothermal, so an energy balance on the reacting fluid, Equation 6.9, is required. There are no shafts or moving boundaries (other than the agitator which is assumed to do negligible work), so the rate of doing work, \(\dot{W}\), is zero. The reacting fluid is a liquid, so the pressure is constant and its time-derivative is zero. Assuming it to be an incompressible ideal mixture means that the volume also is constant, and the time-derivative of the volume is equal to zero. The assignment provides a mass-specific heat capacity, so the sensible heat term can be written in terms of that instead of molar heat capacities. Since there are two reactions, the final summation expands to two terms. After making all those substitutions, both sides of the the reacting fluid energy balance can be divided by \(\rho V \tilde{C}_p\) to put it in the form of a derivative expression.

\[ \cancelto{\rho V \tilde{C}_p}{\left(\sum_i n_i \hat C_{p,i} \right)} \frac{dT}{dt} - \cancelto{0}{V\frac{dP}{dt}} - \cancelto{0}{P \frac{dV}{dt}} = \dot Q - \cancelto{0}{\dot W} - V \cancelto{\left(r_1 \Delta H_1 + r_2 \Delta H_2\right)}{\sum_j \left(r_j \Delta H_j \right)} \]

\[ \rho V \tilde{C}_p \frac{dT}{dt} = \dot Q - V \left(r_1 \Delta H_1 + r_2 \Delta H_2\right) \]

The reactor is cooled with chilled water which gains sensible heat from the reacting fluid, so an energy balance on that exchange fluid is necessary and is given by Equation 6.2.

\[ \rho_{ex} V_{ex} \tilde C_{p,ex}\frac{dT_{ex,out}}{dt} = -\dot Q - \dot m_{ex} \int_{T_{ex,in}}^{T_{ex,out}} \tilde C_{p,ex}dT \]

The exchange fluid heat capacity is constant, making the evaluation of the integral trivial.

\[ \dot m_{ex}\int_{T_{ex,in}}^{T_{ex,out}} \tilde C_{p,ex}dT \Rightarrow \dot m_{ex}\tilde C_{p,ex}\left( T_{ex,out} - T_{ex,in} \right) \]

Dividing both sides of the exchange fluid energy balance by \(\rho_{ex} V_{ex} \tilde C_{p,ex}\) then puts it in the form of a derivative expression. Momentum balances are not used with BSTRS so the full set of design equations includes mole balances on each of the reagents, A, B, X, Y, and Z, an energy balance on the reacting fluid, and an energy balance on the exchange fluid. They consist of seven ODEs containing seven dependent variables, so it isn’t necessary to add an ODE or eliminate a dependent variable.

BSTR Design Equations

Mole balance design equations for A, B, X, Y, and Z are presented in equations (5) through (9). An energy balance on the reacting fluid is given by equation (10), and an energy balance on the exchange fluid is given by equation (11).

\[ \frac{dn_A}{dt} = \left(-r_1 -r_2 \right)V \tag{5} \]

\[ \frac{dn_B}{dt} = -r_1 V \tag{6} \]

\[ \frac{dn_X}{dt} = r_1 V \tag{7} \]

\[ \frac{dn_Y}{dt} = r_1 V \tag{8} \]

\[ \frac{dn_Z}{dt} = r_2 V \tag{9} \]

\[ \frac{dT}{dt} = \frac{\dot Q - V \left(r_1 \Delta H_1 + r_2 \Delta H_2\right)}{\tilde{C}_p \rho V} \tag{10} \]

\[ \frac{dT_{ex,out}}{dt} = \frac{-\dot Q - \dot m_{ex} \tilde C_{p,ex}\left(T_{ex,out} - T_{ex,in}\right)}{\rho_{ex} V_{ex} \tilde C_{p,ex}} \tag{11} \]

Initial values and a stopping ctierion are needed when solving IVODEs or DAEs numerically. The instant the fluids are mixed in the reactor can be defined as \(t=0\). The molar amounts of A, B, X, Y, and Z, the reacting fluid temperature, and the exhhange fluid temperature at that moment are then the initial values. The solutions charged to the reactor did not contain X, Y, or Z, so their initial values are zero. The assignment provides the final value of \(t_f\), which I will use to define the stopping criterion.

The initial values of \(n_A\) and \(n_B\) can be calculated directly from their known initial concentrations and the fluid volume, using the defining equation for concentration in a closed system, Equation 2.7. The initial value of \(T\) is not known and cannot be calculated directly, but a second final value, that of \(n_A\), can, by using its known conversion.

This makes the design equations DAEs because the initial value of one dependent variable is unknown, but two final values are known. As described in Appendix F.5.3, that means I also need to write an implicit ATE relating the known final value of \(n_A\) to the unknown initial temperature. Since the ATE will be solved numerically, I will write it in the form of a residual expression.

Intial Values and Stopping Criterion

Table 9.2: Initial values and stopping criterion for solving design equations, (5) through (11).
Variable Initial Value Stopping Criterion
\(t\) \(0\) \(t_f\)
\(n_A\) \(n_{A,0}\)
\(n_B\) \(n_{B,0}\)
\(n_X\) \(0\)
\(n_Y\) \(0\)
\(n_Z\) \(0\)
\(T\) \(T_0\)
\(T_{ex,out}\) \(T_{ex,out,0}\)

Supporting Equations for Calculating the Initial and Final Values

\[ n_{A,0} = C_{A,0}V \tag{12} \]

\[ n_{B,0} = C_{B,0}V \tag{13} \]

\[ n_{A,f} = n_{A,0} \left( 1 - f_{A,f} \right) \tag{14} \]

Implicit ATE Relating \(T_0\) to \(n_{A,f}\)

\[ T_0 \, : \, n_A\big\vert_{t = t_f} = n_{A,f} \qquad \Rightarrow \qquad \epsilon_{T_0} = n_A\big\vert_{t = t_f} - n_{A,f}\tag{15} \]

To solve the DAEs numerically, I’ll need to write a derivatives function and a residuals function. The residuals function will be passed a guess for \(T_0\) as an argument, and it must return the residual \(\epsilon_{T_0}\). To do so, it must define the initial values and stopping criterion as already described. Then it will call an IVODE solver to solve the BSTR design equations, and use the final value of \(n_A\) that the solver returns to evaluate the residual, equation (15).

The derivatives function will be passed values of the independent and dependent variables, \(t\), \(n_A\), \(n_B\), \(n_X\), \(n_Y\), \(n_Z\), \(T\), and \(T_{ex,out}\), at the start of an integration step as arguments. It must return the derivatives of the dependent variables, calculated using equations (5) through (11). Before it can do that, it will need to calculate any additional unknowns that appear in the derivatives expressions. Upon examination of equations (5) through (11), I see that \(r_1\), \(r_2\), and \(\dot{Q}\) are additional unknowns.

The two rates, \(r_1\) and \(r_2\), can be calculated using the rate expressions, equation (3) and (4). The rate expressions contain four additional variables that will need to be calculated, \(k_1\), \(k_2\), \(C_A\) and \(C_B\). The rate coefficients can be calculated using the Arrhenius expression, Equation 4.8. The concentrations can be calculated using the defining equation for concentration in a closed system, Equation 2.7.

The rate of heat exchange, \(\dot{Q}\), can be calculated using the given heat transfer coefficient, \(U\), and heat transfer area, \(A\), using Equation E.4.

Ancillary Equations for Evaluating the Derivatives

\[ k_1 = k_{0,1}\exp{\left( \frac{-E_1}{RT} \right) } \tag{16} \]

\[ k_2 = k_{0,2}\exp{\left( \frac{-E_2}{RT} \right) } \tag{17} \]

\[ C_A = \frac{n_A}{V} \tag{18} \]

\[ C_B = \frac{n_B}{V} \tag{19} \]

\[ \dot{Q} = UA\left(T_{ex,out} - T \right) \tag{20} \]

At this point the DAEs, equations (5) through (11) and (15), can be solved numerically. Doing so will yield \(T_0\) along with corresponding sets of values of \(t\), \(n_A\), \(n_B\), \(n_X\), \(n_Y\), \(n_Z\), \(T\), and \(T_{ex,out}\) that span the range from their initial to their final values. \(T_0\) and the final values of \(T\), and \(T_{ex,out}\) are quantities of interest. The other quantity of interest, the final selectivity for X relative to Z, can be calculated using the defining equation for selectivity, Equation 3.17.

9.6.2.3 Calculating the Quantities of Interest

With the information given above, the DAEs, equations (5) through (11) and (15), can be solved numerically, yielding \(T_0\) along with corresponding sets of values of \(t\), \(n_A\), \(n_B\), \(n_X\), \(n_Y\), \(n_Z\), \(T\), and \(T_{ex,out}\) that span the range from their initial to their final values. \(T_0\) is one of the quantities of interest. The final temperature, exchange fluid temperature and selectivity can be calculated using equations (21), (22), and (23).

\[ T_f = T\big\vert_{t=t_f} \tag{21} \]

\[ T_{ex,out,f} = T_{ex,out}\big\vert_{t=t_f} \tag{22} \]

\[ S_{X/Z,f} = \frac{n_X \Big\vert_{t=t_{rxn}}}{n_Z \Big\vert_{t=t_{rxn}}} \tag{23} \]

9.6.2.4 Formulation of the Calculations

To complete the assignment the BSTR design equations are solved first, and then the quantities of interest are calculated. The BSTR design equations are DAEs, and they can be solved numerically as described in Appendix F.5.3 by writing a residuals function, a derivatives function, and a BSTR function. A deliverables function can be added to calculate the quantities of interest.

  • The residuals function must receive a guess for the ATE variable, \(T_0\), and it must return the residual, \(\epsilon_{T_0}\). To do so it must define the initial values and stopping criterion for the ODEs, pass them to an IVODE solver along with the derivatives function, and use the results it returns to evaluate the residual.
  • The derivatives function must receive the independent and dependent variables at the start of an integration step as arguments, and it must return the derivatives of the dependent variables. To do so it must calculate the rate and the rate of heat transfer before it can evaluate the derivatives.
  • For this assignment, the BSTR function must receive an initial guess for the ATE variable and it must return the initial temperature, \(T_0\), and sets of corresponding values of the independent and dependent variables spanning the interval from their initial to their final values. To do so it must pass the guess for \(T_0\) to an ATE solver along with the residuals function, use the value of \(T_0\) returned by the solver to define the initial values and stopping criterion for the ODEs, and then pass them to an IVODE solver along with the derivatives function.
  • The deliverables function does not need any arguments, and it doesn’t return anything. It defines an initial guess for the ATE variable, calls the BSTR function and uses the results it returns to calculate the quantities of interest.

Specifications can be written for each of these functions in preparation for writing them. The function specifications given here assume that the known constants are available wherever they are needed.

Residuals Function Specifications

Argument: \(T_0\)

Returns: \(\epsilon_{T_0}\)

Algorithm:

\[ t_{initial} = 0 \tag{24} \]

\[ n_{A,initial} = n_{A,0} = C_{A,0}V \tag{12} \]

\[ n_{B,initial} = n_{B,0} = C_{B,0}V \tag{13} \]

\[ \begin{matrix} n_{X,initial} = 0,\, n_{Y,initial} = 0,\, n_{Z,initial} = 0 \\ T_{initial} = T_0,\, T_{ex,out,initial} = T_{ex,out,0},\, t_{final} = t_f \end{matrix} \tag{25} \]

 

\[ \begin{matrix} t_{initial}, n_{A,initial}, n_{B,initial}, n_{X,initial}, n_{Y,initial}, \\ n_{Z,initial}, T_{initial}, n_{A,final}, T_{ex,out,initial}, \\ \text{derivatives function}\\ \Downarrow \\ \text{IVODE solver} \\ \Downarrow \\ \underline{t}, \underline{n}_A, \underline{n}_B, \underline{n}_Y, \underline{n}_Z, \underline{T} ,\underline{T}_{ex,out} \end{matrix} \tag{26} \]

 

\[ n_{A,f} = n_{A,0} \left( 1 - f_{A,f} \right) \tag{14} \]

\[ \epsilon_{T_0} = n_A\big\vert_{t = t_f} - n_{A,f}\tag{15} \]

Derivatives Function Specifications

Arguments: \(t\), \(n_A\), \(n_B\), \(n_X\), \(n_Y\), \(n_Z\), \(T\), and \(T_{ex,out}\)

Returns: \(\frac{dn_A}{dt}\), \(\frac{dn_B}{dt}\), \(\frac{dn_X}{dt}\), \(\frac{dn_Y}{dt}\), \(\frac{dn_Z}{dt}\), \(\frac{dT}{dt}\), and \(\frac{dT_{ex,out}}{dt}\)

Algorithm:

\[ k_1 = k_{0,1}\exp{\left( \frac{-E_1}{RT} \right) } \tag{16} \]

\[ k_2 = k_{0,2}\exp{\left( \frac{-E_2}{RT} \right) } \tag{17} \]

\[ C_A = \frac{n_A}{V} \tag{18} \]

\[ C_B = \frac{n_B}{V} \tag{19} \]

\[ \dot{Q} = UA\left(T_{ex,out} - T \right) \tag{20} \]

\[ r_1 = k_1C_AC_B \tag{3} \]

\[ r_2 = k_2C_A \tag{4} \]

\[ \frac{dn_A}{dt} = \left(-r_1 -r_2 \right)V \tag{5} \]

\[ \frac{dn_B}{dt} = -r_1 V \tag{6} \]

\[ \frac{dn_X}{dt} = r_1 V \tag{7} \]

\[ \frac{dn_Y}{dt} = r_1 V \tag{8} \]

\[ \frac{dn_Z}{dt} = r_2 V \tag{9} \]

\[ \frac{dT}{dt} = \frac{\dot Q - V \left(r_1 \Delta H_1 + r_2 \Delta H_2\right)}{\tilde{C}_p \rho V} \tag{10} \]

\[ \frac{dT_{ex,out}}{dt} = \frac{-\dot Q - \dot m_{ex} \tilde C_{p,ex}\left(T_{ex,out} - T_{ex,in}\right)}{\rho_{ex} V_{ex} \tilde C_{p,ex}} \tag{11} \]

BSTR Function Specifications

Arguments: \(T_{0,guess}\)

Returns: \(T_0\), \(\underline{t}\), \(\underline{n}_A\), \(\underline{n}_B\), \(\underline{n}_Y\), \(\underline{n}_Z\), \(\underline{T}\) ,\(\underline{T}_{ex,out}\)

Algorithm:

\[ \begin{matrix} T_{0,guess},\, \text{residuals function}\\ \Downarrow \\ \text{ATE solver} \\ \Downarrow \\ T_0 \end{matrix}\tag{27} \]

 

\[ t_{initial} = 0 \tag{24} \]

\[ n_{A,initial} = n_{A,0} = C_{A,0}V \tag{12} \]

\[ n_{B,initial} = n_{B,0} = C_{B,0}V \tag{13} \]

\[ \begin{matrix} n_{X,initial} = 0,\, n_{Y,initial} = 0,\, n_{Z,initial} = 0 \\ T_{initial} = T_0,\, T_{ex,out,initial} = T_{ex,out,0},\, t_{final} = t_f \end{matrix} \tag{25} \]

 

\[ \begin{matrix} t_{initial}, n_{A,initial}, n_{B,initial}, n_{X,initial}, n_{Y,initial}, \\ n_{Z,initial}, T_{initial}, n_{A,final}, T_{ex,out,initial}, \\ \text{derivatives function}\\ \Downarrow \\ \text{IVODE solver} \\ \Downarrow \\ \underline{t}, \underline{n}_A, \underline{n}_B, \underline{n}_Y, \underline{n}_Z, \underline{T} ,\underline{T}_{ex,out} \end{matrix} \tag{26} \]

The deliverables function must define an initial guess for \(T_0\) to pass as an argument to the BSTR function. The only temperatures I know are the inlet exchange fluid temperature and its initial outlet temperature, and they are both equal to 40 °C. So I’ll guess that the initial reacting fluid temperature is also 40 °C, recognizing that if the ATE solver fails to converge, I’ll need to repeat the calculations using a different initial guess.

Deliverables Function Specifications

Arguments: none

Returns: nothing

Algorithm:

\[ T_{0,guess} = \text{40 °C} \tag{28} \]

 

\[ \begin{matrix} T_{0,guess}\\ \Downarrow \\ \text{BSTR function} \\ \Downarrow \\ T_0, \underline{t}, \underline{n}_A, \underline{n}_B, \underline{n}_Y, \underline{n}_Z, \underline{T} ,\underline{T}_{ex,out} \end{matrix} \tag{29} \]

 

\[ T_f = T\big\vert_{t=t_f} \tag{21} \]

\[ T_{ex,out,f} = T_{ex,out}\big\vert_{t=t_f} \tag{22} \]

\[ S_{X/Z,f} = \frac{n_X \Big\vert_{t=t_{rxn}}}{n_Z \Big\vert_{t=t_{rxn}}} \tag{23} \]

 

\[ \begin{matrix} T_{0}, T_f, T_{ex,out,f}, S_{X/Z,f}\\ \Downarrow \\ \text{functions to display and save} \\ \end{matrix} \tag{30} \]

9.6.2.5 Implementation and Execution of the Calculations

The residuals, derivatives, BSTR, and deliverables functions were written as specified above using an appropriate mathematics software package. The calculations were performed by calling the deliverables function. Upon doing so, the ATE solver converged, so it was not necessary to repeat the calculations using a different initial guess for the ATE variable, \(T_0\).

9.6.2.6 Results and Discussion

The initial temperature must be 65 °C for the conversion to equal 45% after 30 min of reaction. At that time, the temperature of the reacting fluid will be 92.4 °C, the exit temperature of the exchange fluid will be 68.2 °C, and the selectivity will equal 4.21 mol X per mol Z.

Solving the reactor design equations using an initial temperature of 65 °C also yields the final molar amount of \(n_A\). One way to check the results is to use that value to calculate the conversion and make sure that it is equal to 45%. Indeed, when the design equations were solved using an initial temperature of 65 °C, the resulting conversion of A was 45%.

When two reactions are taking place, it is quite common for the product of one reaction to be desired and the product of the other reaction to be undesired. That raises the question, “what can a reaction engineer do to change the selectivity?” There isn’t a single, universal answer to that question. However, for parallel reactions it is often useful to write an expression for the instantaneous selectivity using Equation 4.5. In this system, the instantaneous selectivity for X over Z is given by the net rate of generation of X divided by the net rate of generation of Z as shown in equation (31).

\[ \begin{align} S_{X/Z,inst} &= \frac{\displaystyle \sum _j r_{X,j}}{\displaystyle \sum _j r_{Z,j}} = \frac{r_1}{r_2} = \frac{k_1C_AC_B}{k_2C_A} = \frac{k_{0,1}\exp{\left(\frac{-E_1}{RT}\right)}}{k_{0,2}\exp{\left(\frac{-E_2}{RT}\right)}}C_B \\ S_{X/Z,inst} &= \frac{k_{0,1}}{k_{0,2}}\exp{\left(\frac{E_2 - E_1}{RT}\right)}C_B \end{align} \tag{31} \]

The final selectivity will equal the integral of the instantaneous selectivity from \(t_0\) to \(t_f\). Noting that \(E_2\) is greater than \(E_1\), equation (31) suggests two ways to increase the instantantous selectivity for X over Z. The first is to maintain a high concentration of B. The second is to decrease the mean temperature during the reaction. In a batch reactor, the only way to increase the average concentration of B is to increase the initial amount. The average temperature could be decreased by using a lower initial temperature or by changing the exchange fluid flow rate or inlet temperature. One point to bear in mind is that decreasing the mean temperature will decrease both reaction rates, not just the undesired rate. As a consequence, to keep the conversion the same, a longer reaction time would be necessary.

To test these predictions, a simulation was performed where the initial temperature was reduced to 55 °C, keeping the conversion the same at 45%. The results showed that the selectivity will increase to 5.51 mol X per mol Z, but the reaction time will increase to 87.3 min. Given economic data and other operating constraints, an engineer could optimize the performance of this reactor through simulation of this kind.

Note

This example specified two final values: the reaction time and the fractional conversion. In the solution presented above, the reaction time was used as the stopping criterion when solving the IVODEs and the fractional conversion was used in the implicit equation for \(T_0\).

The roles of the reaction time and conversion could have been reversed. The specified conversion could have been used to calculate the corresponding final moles of A, which then could have been used as the stopping criterion. Then the reaction time would have been used in the implicit equation for \(T_0\). The results would be the same.

9.6.3 Maximum Intermediate Yield in Series Reactions

A BSTR will be charged with 1 atm of A and 2 atm of B at 25 °C. The reactor volume is 2 L, and it has a jacket with an area of 600 cm2 and an overall heat transfer coefficient of 0.6 cal cm-2 min-1 K-1. The temperature of the coolant in the perfectly-mixed jacket is maintained at 30°C by a temperature controller. Gas phase reactions (1) and (2) take place within the reactor; the corresponding rate expressions are given in equations (3) and (4). The rate coefficients display Arrhenius temperature dependence with pre-exponental factors of 3.34 x 109 and 1.47 x 1010 mol cm-3 min-1 atm-2 for reactions (1) and (2), respectively, and activation energies of 20.5 and 21.8 kcal mol-1. The heat of reaction (1) is constant and equal to -6,300 cal mol-1; that of reaction (2) is constant and equal to -6,900 cal mol-1. The heat capacities of A, B, D, Z and U are constant and equal to 7.4, 8.6, 10.7, 5.2 and 10.3 cal mol-1 K-1, respectively.

What reaction time will maximize the yield (moles of D per initial mole of A), and what will the yield and the conversion of A equal at that reaction time?

\[ A + B \rightarrow D + Z \tag{1} \]

\[ D + B \rightarrow U + Z \tag{2} \]

\[ r_1 = k_{0,1} \exp{\left(\frac{-E_1}{RT}\right)}P_AP_B \tag{3} \]

\[ r_2 = k_{0,2} \exp{\left(\frac{-E_2}{RT}\right)}P_DP_B \tag{4} \]


This assignment describes a single BSTR and no other equipment. It describes the reactor, its operation, and the reactions taking place, including the rate expressions. It asks me to find a reaction time to maximize yield. From these characteristics, I know that this is an isolated reactor modeling assignment. I’ll begin by summarizing the assignment. I’ll use appropriate variable symbols to represent the given quantities and the quantities of interest. A cooled, batch reactor is being used, so I’ll use a subscripted “0” to denote initial values and a subscripted “f” to denote final values. I’ll let \(t_{opt}\) denote the reaction time that maximizes the yield. The exchange fluid temperature is constant and known, but no information is provided about the inlet temperature or flow rate of the heat exchange fluid, so there’s no need for subscripts denoting “in” and “out.”

9.6.3.1 Assignment Summary

Reactions:

\[ A + B \rightarrow D + Z \tag{1} \]

\[ D + B \rightarrow U + Z \tag{2} \]

Rate Expressions:

\[ r_1 = k_{0,1} \exp{\left(\frac{-E_1}{RT}\right)}P_AP_B \tag{3} \]

\[ r_2 = k_{0,2} \exp{\left(\frac{-E_2}{RT}\right)}P_DP_B \tag{4} \]

Reactor System: BSTR with heat transfer

Reactor Schematic:

Figure 9.8: A closed (batch), non-adiabatic stirred tank reactor with a constant exchange fluid temperature.

Quantities of Interest: \(t_{opt} = \underset{t}{\arg\max} \left(Y_{D/A}\right)\), \(Y_{D/A}\Big\vert_{t_{opt}}\), and \(f_A\Big\vert_{t_{opt}}\).

Given and Known Constants: \(P_{A,0}\) = 1 atm, \(P_{B,0}\) = 2 atm, \(T_0\) = 25 °C, \(V\) = 2 L, \(A\) = 600 cm2, \(U\) = 0.6 cal cm-2 min-1 K-1, \(T_{ex}\) = 30 °C, \(k_{0,1}\) = 3.34 x 109 mol cm-3 min-1 atm-2, \(k_{0,2}\) = 1.47 x 1010 mol cm-3 min-1 atm-2, \(E_1\) = 20.5 kcal mol-1, \(E_2\) = 21.8 kcal mol-1, \(\Delta H_1\) = -6,300 cal mol-1, \(\Delta H_2\) = -6,900 cal mol-1, \(\hat C_{p,A}\) = 7.4 cal mol-1 K-1, \(\hat C_{p,B}\) = 8.6 cal mol-1 K-1, \(\hat C_{p,D}\) = 10.7 cal mol-1 K-1, \(\hat C_{p,Z}\) = 5.2 cal mol-1 K-1, \(\hat C_{p,U}\) = 10.3 cal mol-1 K-1, and \(R\) = 1.987 cal mol-1 K-1 = 82.057 cm3 atm mol-1 K-1.

9.6.3.2 Generation of the Model Equations

Mole balances are always included in the reactor design equations. The general BSTR mole balance is given in Equation 6.8. I will use this equation to write mole balances for every reagent in the system. That is, I’ll write the mole balance for \(i\) = A, B, D, Z, and U. There are two reactions taking place, so the sum will expand to two terms in each mole balance.

\[ \frac{dn_i}{dt} = V \sum_j \nu_{i,j}r_j = \left( \nu_{i,1}r_1 + \nu_{i,2}r_2 \right)V \]

This reactor is not isothermal, so I will need to include a BSTR reacting fluid energy balance among the reactor design equations. The BSTR energy balance is given in Equation 6.9. This reactor has rigid walls and no moving boundaries, so it’s volume is constant. If the volume is constant, the time-derivative of the volume is equal to zero. Being a gas phase system with a constant volume, the pressure is expected to change because the temperature will change. Therefore I cannot set the time-derivative of the pressure equal to zero. With no moving boundaries, the only work is that of the agitator, which is assumed to be negligible. The summation over \(i\) includes every reagent, and as above, the summation over \(j\) expands to two terms.

\[ \left(\sum_i n_i \hat C_{p,i} \right) \frac{dT}{dt} - V\frac{dP}{dt} - P \cancelto{0}{\frac{dV}{dt}} = \dot Q - \cancelto{0}{\dot W} - V \cancelto{\left( r_1 \Delta H_1 + r_2 \Delta H_2\right)}{\sum_j \left(r_j \Delta H_j \right)} \]

\[ \left( n_A \hat C_{p,A} + n_B \hat C_{p,B} + n_D \hat C_{p,D} + n_Z \hat C_{p,Z} + n_U \hat C_{p,U} \right) \frac{dT}{dt} - V\frac{dP}{dt} = \dot Q - \left( r_1 \Delta H_1 + r_2 \Delta H_2\right)V \]

An exchange fluid is present in this reactor, but the assignment states that its temperature is constant at 30 °C. Knowing the exchange fluid temperature, the other reactor design equations can be solved independently of an energy balance on the exchange fluid, and because the assignment does not ask anything about the exchange fluid, I do not need to include an energy balance on the exchange fluid among the reactor design equations.

BSTR Design Equations

Mole balances of each of the reagents are presented in equations (5) through (9), and an energy balance on the reacting fluid is presented in equation (10).

\[ \frac{dn_A}{dt} = -r_1V \tag{5} \]

\[ \frac{dn_B}{dt} = \left( -r_1 - r_2 \right)V \tag{6} \]

\[ \frac{dn_D}{dt} = \left( r_1 - r_2 \right)V \tag{7} \]

\[ \frac{dn_Z}{dt} = \left( r_1 + r_2 \right)V \tag{8} \]

\[ \frac{dn_U}{dt} = r_2V \tag{9} \]

\[ \begin{aligned} - V\frac{dP}{dt} + &\left( n_A \hat C_{p,A} + n_B \hat C_{p,B} + n_D \hat C_{p,D} + n_Z \hat C_{p,Z} + n_U \hat C_{p,U} \right) \frac{dT}{dt} \\&= \dot Q - \left( r_1 \Delta H_1 + r_2 \Delta H_2\right)V \end{aligned} \tag{10} \]

At this point I have 6 initial value ordinary differential equations that contain 7 dependent variables (\(n_A\), \(n_B\), \(n_D\), \(n_Z\), \(n_U\), \(P\), and \(T\)). I either need to eliminate a dependent variable or add an ODE before I can solve the reactor design equations.

Noting that V is constant, I’m going to take the time-derivative of the ideal gas law and use it as an additional ODE.

\[ PV - \left(n_A + n_B + n_D + n_Z + n_U \right)RT =0 \]

\[ \begin{aligned} V\frac{dP}{dt} &- RT\left(\frac{dn_A}{dt} + \frac{n_B}{dt} + \frac{n_D}{dt} + \frac{n_Z}{dt} + \frac{n_U}{dt} \right) \\&- R\left(n_A + n_B + n_D + n_Z + n_U \right)\frac{dT}{dt} = 0 \end{aligned} \]

Taking the derivative of the ideal gas law yields equation (11).

\[ \begin{aligned} - &RT\left(\frac{dn_A}{dt} + \frac{n_B}{dt} + \frac{n_D}{dt} + \frac{n_Z}{dt} + \frac{n_U}{dt} \right) + V\frac{dP}{dt} \\&- R\left(n_A + n_B + n_D + n_Z + n_U \right)\frac{dT}{dt} = 0 \end{aligned}\tag{11} \]

The reactor design equations are either IVODEs or DAEs, so initial values and a stopping criterion are needed in order to solve them. I can define \(t=0\) to be instant that the A and B are added to the reactor. The other initial values are then the molar amounts of each reagent, the temperature, and the pressure at that time. The assignment states that the BSTR is charged with only A and B, so the initial molar amounts of D, Z, and U are zero.

The initial molar amounts of A and B can be calculated directly using their known initial partial pressures and the ideal gas law. The initial total pressure is simply the sum of the initial partial pressures.

I need to find the value of \(t_f\) that maximizes the yield of D. To do that, I’ll choose a large value for the final time, \(t_f\). Upon doing so, all of the initial values and one final value will be known or can be calculated directly, so the BSTR design equations in this assignment are IVODEs. As long as the final value I choose is large enough, the yield will reach its maximum before the final time is reached. The time when it reaches its maximum is the optimum time requested in the assignment narrative.

Initial Values and Stopping Criterion

Table 9.3: Initial values and stopping criterion for solving the design equations, (5) through (11).
Variable Initial Value Stopping Criterion
\(t\) \(0\) \(t_f\)
\(n_A\) \(n_{A,0}\)
\(n_B\) \(n_{B,0}\)
\(n_D\) \(0\)
\(n_Z\) \(0\)
\(n_U\) \(0\)
\(T\) \(T_0\)
\(P\) \(P_0\)

Supporting Equations for Calculating the Initial and Final Values

\[ n_{A,0} = \frac{P_{A,0}V}{RT_0} \tag{12} \]

\[ n_{B,0} = \frac{P_{B,0}V}{RT_0} \tag{13} \]

\[ P_0 = P_{A,0} + P_{B,0} \tag{14} \]

I now have 7 IVODEs with 7 dependent variables. In order to solve the IVODEs numerically I’ll need to calculate the values of the derivatives at the start of each integration step. To begin, I need to write the IVODEs in the form of derivative expressions. I’ll use the procedure described in Appendix F.4.2. It involves writing the IVODEs in the form of a matrix equation.

\[ \boldsymbol{M}\frac{d}{dt} \begin{bmatrix} n_A \\ n_B \\ n_D \\ n_Z \\ n_U \\ T \\ P \end{bmatrix} = \underline{g} \]

In other words, I’ll need to define the elements of the mass matrix, \(\boldsymbol{M}\), and the vector, \(\underline{g}\). If anything other than known constants, the independent variable, and the dependent variables appear in the elements of \(\boldsymbol{M}\) and \(\underline{g}\) I’ll need to calculate them before I can evaluate the derivatives. Looking at the design equations, I see that I will need to calculate \(r_1\), \(r_2\), and \(\dot{Q}\).

The rate expressions, equations (3) and (4), can be used to calculate the rates, but those equations introduce the additional variables, \(P_A\), \(P_B\), and \(P_D\). Those partial pressures can be calculated using the ideal gas law.

The rate of heat transfer, \(\dot{Q}\), can be calculated using the known heat transfer coefficient and heat transfer area.

Ancillary Equations for Evaluating the Derivatives

The IVODEs, equations (1) through (7) can be written as the matrix equation, (15), with the mass matrix defined as shown in equation (16) and with the vector, \(\underline{g}\), defined as shown in equation (17). The IVODEs can then be written as derivative expressions, equation (18).

\[ \boldsymbol{M}\frac{d}{dt} \begin{bmatrix} n_A \\ n_B \\ n_D \\ n_Z \\ n_U \\ T \\ P \end{bmatrix} = \underline{g}\tag{15} \]

\[ \boldsymbol{M} = \begin{bmatrix} 1 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & \displaystyle\sum_i\left(n_i \hat{C}_{p,i} \right) & -V \\ RT & RT & RT & RT & RT & R\displaystyle\sum_i n_i & -V \end{bmatrix} \tag{16} \]

\[ \underline{g} = \begin{bmatrix} -r_1V \\ \left(-r_1 - r_2\right)V \\ \left(r_1 - r_2\right)V \\ \left(r_1 + r_2\right)V \\ r_2V \\ \dot{Q} -\left(r_1\Delta H_1 + r_2 \Delta H_2\right)V \\ 0 \end{bmatrix}\tag{17} \]

\[ \begin{bmatrix} \frac{dn_A}{dt} \\ \frac{n_B}{dt} \\ \frac{n_D}{dt} \\ \frac{n_Z}{dt} \\ \frac{dn_U}{dt} \\ \frac{dT}{dt} \\ \frac{dP}{dt} \end{bmatrix} = \boldsymbol{M}^{-1} \underline{g} \tag{18} \]

The rates, \(r_1\) and \(r_2\), can be calculated using equations (3) and (4). The partial pressures appearing in the rate expressions can be calculated using equations (19) through (21), and the rate of heat exchange can be calculted using equation (22)

\[ P_A = \frac{n_ART}{V} \tag{19} \]

\[ P_B = \frac{n_BRT}{V} \tag{20} \]

\[ P_D = \frac{n_DRT}{V} \tag{21} \]

\[ \dot Q = UA\left( T_{ex} - T \right) \tag{22} \]

At this point the BSTR design equations can be solved for sets of values of \(\underline{t}\), \(\underline{n}_A\), \(\underline{n}_B\), \(\underline{n}_D\), \(\underline{n}_Z\), \(\underline{n}_U\), \(\underline{T}\), and \(\underline{P}\) spanning the range from \(t_0\) to the chosen final time, \(t_f\). At each time the yield can be calculated using Equation 3.13, and the conversion can be calculated using the defining equation for conversion in a closed system, Equation 3.4.

The resulting set of yields can then be examined to find the maximum value and the time, \(t_{opt}\), when it occurred. If the yield does not pass through a maximum, the calculations must be repeated using a larger chosen final time.

9.6.3.3 Calculating the Quantities of Interest

The BSTR design equations, (5) through (11), can be solved numerically to obtain corresponding sets of values of \(t\), \(n_A\), \(n_B\), \(n_D\), \(n_Z\), \(n_U\), \(P\), and \(T\) that span the range from \(t_0\) to \(t_f\). At each time, \(t_i\), the corresponding yield of D, \(Y_{D/A,i}\), and conversion of A, \(f_{A,i}\), can be calculated using equations (23) and (24). The maximum yield, the time at which it occurs, and the corresponding conversion of A can then be identified, equations (25) through (28). If the yield does not pass through a maximum, the calculations must be repeated using a larger chosen final time.

\[ Y_{D/A,i} = \frac{n_{D,i}}{n_{A,0}} \tag{23} \]

\[ f_{A,i} = \frac{n_{A,0} - n_{A,i}}{n_{A,0}} \tag{24} \]

\[ i_{opt} = \underset{i}{\arg\max} \left(Y_{D/A,i}\right) \tag{25} \]

\[ Y_{D/A}\big\vert_{t_{opt}} = Y_{D/A,i_{opt}} \tag{26} \]

\[ t_{opt} = t_{f,i_{opt}}\tag{27} \]

\[ f_A\big\vert_{t_{opt}} = f_{A,i_{opt}} \tag{28} \]

9.6.3.4 Formulation of the Calculations

To complete the assignment, the BSTR design equations are solved first, and then the quantities of interest are calculated. The BSTR design equations are IVODEs, and they can be solved numerically as described in Appendix F.5.2 by writing a derivatives function and a BSTR function. A deliverables function can be added to calculate the quantities of interest.

  • The derivatives function must be passed the independent and dependent variables at the start of an integration step, and it must return the derivatives of the dependent variables. To do so the algorithm must calculate the rates, \(r_1\) and \(r_2\), and the rate of heat transfer \(\dot{Q}\), form the mass matrix, \(\boldsymbol{M}\) and the vector \(g\), and evaluate the derivatives.

  • The BSTR function will be passed a final time as an argument, and return corresponding sets of the independent and dependent variables spanning the range from their initial to their final values. The algorithm will simply define the initial values and stopping criterion and then call an IVODE solver.

  • The deliverables function will take no arguments and return nothing. It will define a large final time, call the BSTR function, use the results to calculate the quantities of interest, display, and save them to a file.

Specifications can be written for each of these functions in preparation for writing them. The function specifications given here assume that the known constants are available wherever they are needed.

Derivatives Function Specifications

Arguments: \(t\), \(n_A\), \(n_B\), \(n_D\), \(n_Z\), \(n_U\), \(P\), and \(T\)

Returns: \(\frac{dn_A}{dt}\), \(\frac{dn_B}{dt}\), \(\frac{dn_D}{dt}\), \(\frac{dn_Z}{dt}\), \(\frac{dn_U}{dt}\), \(\frac{dP}{dt}\), and \(\frac{dT}{dt}\)

Algorithm:

\[ P_A = \frac{n_ART}{V} \tag{19} \]

\[ P_B = \frac{n_BRT}{V} \tag{20} \]

\[ P_D = \frac{n_DRT}{V} \tag{21} \]

\[ r_1 = k_{0,1} \exp{\left(\frac{-E_1}{RT}\right)}P_AP_B \tag{3} \]

\[ r_2 = k_{0,2} \exp{\left(\frac{-E_2}{RT}\right)}P_DP_B \tag{4} \]

\[ \dot Q = UA\left( T_{ex} - T \right) \tag{22} \]

\[ \boldsymbol{M} = \begin{bmatrix} 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & \displaystyle\sum_i\left(n_i \hat{C}_{p,i} \right) & -V \\ RT & RT & RT & RT & RT & R\displaystyle\sum_i n_i & -V \end{bmatrix} \tag{16} \]

\[ \underline{g} = \begin{bmatrix} -r_1V \\ \left(-r_1 - r_2\right)V \\ \left(r_1 - r_2\right)V \\ \left(r_1 + r_2\right)V \\ r_2V \\ \dot{Q} -\left(r_1\Delta H_1 + r_2 \Delta H_2\right)V \\ 0 \end{bmatrix}\tag{17} \]

\[ \begin{bmatrix} \frac{dn_A}{dt} \\ \frac{n_B}{dt} \\ \frac{n_D}{dt} \\ \frac{n_Z}{dt} \\ \frac{dn_U}{dt} \\ \frac{T}{dt} \\ \frac{P}{dt} \end{bmatrix} = \boldsymbol{M}^{-1} \underline{g} \tag{18} \]

BSTR Function Specifications

Arguments: \(t_f\)

Returns: \(\underline{t}\), \(\underline{n}_A\), \(\underline{n}_B\), \(\underline{n}_D\), \(\underline{n}_Z\), \(\underline{n}_U\), \(\underline{T}\),and \(\underline{P}\)

Algorithm:

\[ t_{initial} = 0 \tag{33} \]

\[ n_{A,initial} = n_{A,0} = \frac{P_{A,0}V}{RT} \tag{12} \]

\[ n_{B,initial} = n_{B,0} = \frac{P_{B,0}V}{RT} \tag{13} \]

\[ n_{D,initial} = 0,\, n_{Z,initial} = 0,\, n_{U,initial} = 0,\, T_{initial} = T_0 \tag{34} \]

\[ P_{initial} = P_0 = P_{A,0} + P_{B,0} \tag{14} \]

\[ t_{final} = t_f \tag{35} \]

\[ \begin{matrix}t_{initial}, n_{A,initial}, n_{B,initial}, n_{D,initial}, \\ n_{Z,initial}, n_{U,initial}, n_{A,initial}, T_{initial},\\ P_{initial}, t_{final}, \text{ derivatives function}\\ \Downarrow\\ \text{IVODE solver}\\ \Downarrow\\ \underline{t}, \underline{n}_A, \underline{n}_B, \underline{n}_D, \underline{n}_Z, \underline{n}_U, \underline{T},\underline{P}\\ \end{matrix} \tag{36} \]

 

The deliverables function must define a final time that is greater than the time when the yield reaches its maximum, but I don’t know what that time is. So here I’ll set the final time to be 20 minutes. If it takes longer for the yield to reach its maximum, I’ll need to increase the final time and repeat the calculations.

Deliverables Function Specifications

Arguments: none

Returns: nothing

Algorithm:

\[ t_f = 20 \tag{37} \]

 

\[ \begin{matrix}t_f\\ \Downarrow\\ \text{BSTR function}\\ \Downarrow\\ \underline{t}, \underline{n}_A, \underline{n}_B, \underline{n}_D, \underline{n}_Z, \underline{n}_U, \underline{T},\underline{P} \end{matrix} \tag{38} \]

 

\[ \forall\, i \quad Y_{D/A,i} = \frac{n_{D,i}}{n_{A,0}} \tag{23} \]

\[ \forall\, i \quad f_{A,i} = \frac{n_{A,0} - n_{A,i}}{n_{A,0}} \tag{24} \]

\[ i_{opt} = \underset{i}{\arg\max} \left(Y_{D/A,i}\right) \tag{25} \]

\[ Y_{D/A}\big\vert_{t_{opt}} = Y_{D/A,i_{opt}} \tag{26} \]

\[ t_{opt} = t_{f,i_{opt}}\tag{27} \]

\[ f_A\big\vert_{t_{opt}} = f_{A,i_{opt}} \tag{28} \]

 

\[ \begin{matrix} t_{opt},\, Y_{D/A}\big\vert_{t_{opt}},\, f_A\big\vert_{t_{opt}}\\ \Downarrow\\ \text{functions to display and save} \end{matrix} \tag{33} \]

9.6.3.5 Implementation and Execution of the Calculations

The derivatives, BSTR, and deliverables functions were written as specified above using an appropriate mathematics software package. Upon attempting to perform the calculations by calling the deliverables function, the code ran for a long time, and the results were not realistic. This suggested that the IVODEs might be stiff, so the calculations were repeated after setting the code to use a stiff IVODE solver. The code then ran successfully. The yield of D did pass through a maximum, so it was not necessary to change the chosen value of \(t_f\) and repeat the calculations.

9.6.3.6 Results and Discussion

Figure 9.9 shows the variation of the yield as a function of the reaction time. The maximum yield, 0.495 mol D per initial mol A, occurs at a reaction time of 7.1 min. The conversion of A at that reaction time is 74.7 %.

Figure 9.9: Yield of D as the reaction time varies.

Reactions (1) and (2) constitute a series-parallel reaction network. The desired product, D, is an intermediate product. It is produced in reaction (1) and consumed in reaction (2). The yield of an intermediate product is expected to pass through a maximum as the reaction time increases as seen in Figure 9.9. This can be predicted qualitatively as long as the effect of temperature is comparable for the two reactions.

At the start of the reaction, the concentrations of A and B will be large and those of D, U, and Z will be zero. The rate of reaction (1) will be positive while that for reaction (2) will be zero. Consequently, during a very brief interval at the start of the reactive process, the concentrations of A and B will decrease and the concentrations of D and Z will increase.

During the next brief interval in time, the rate of reaction (1) will be smaller because the concentrations of A and B will be smaller, and the rate of reaction (2) no longer will equal zero because a small amount D will be present leading to a positive rate. The concentrations of A and B will again decrease and the concentration of Z will again increase. It can be expected that the rate of reaction (1) will be greater than the rate of reaction (2) because only a small concentration of D will be present, so the concentration of D, and its yield, will again increase, as can be seen in Figure 9.9.

For as long as the rate of reaction (2) remains smaller than the rate of reaction (1), the amount of D will be increasing, but its rate of increase will get smaller and smaller. The reason for this is that the concentrations of A and B will be decreasing, and hence the rate of generation of D by reaction (1) will be decreasing. At the same time, the rate of consumption of D by reaction (2) will be increasing as the concentration of D builds. This can be seen in Figure 9.9 where the curvature of the yield plot is initially concave downward,

Eventually, the rate of reaction (2) will become equal to the rate of reaction (1). For that instant in time, the yield of D will not be changing. This corresponds to the maximum in Figure 9.9. At times beyond the maximum the rate of reaction (1) is less than the rate of reaction (2), so the concentration of D decreases. As the concentrations of both reactants, A and D, continue to decrease, both reactions (1) and (2) become slower and slower. This can be seen at larger reaction times in Figure 9.9 where the curve is concave upward. Eventually, at very large reaction time, the rates will equal zero and all of the D will have been consumed.

In this example series-parallel reactions are occurring where the desired product is produced in one reaction and consumed by another. This assignment shows that the reaction time has a significant effect upon the amount of the desired product. Other process parameters such as the initial composition or temperature will also affect the results, but for series and series-parallel reactions their effect on the selectivity is often weaker. Temperature can be an effective parameter for adjusting the selectivity of multiple reactions if the activation energies of the reactions differ significantly. That was the case in Example 9.6.2, and indeed, the discussion following that example showed that changing the temperature did have a significant effect.

Notes
  1. In the solution above, the mass matrix was used to write the IVODEs in the form of seven derivatives expressions. An equivalent approach is to generate the 7 derivative expressions by algebraic manipulation of equations (5) through (11).

  2. Yet another approach is to solve equation (11) for \(\frac{dP}{dt}\) and substitute the result into equation (10). This reduces the IVODEs to a set of six equations containing six dependent variables: \(n_A\), \(n_B\), \(n_D\), \(n_Z\), \(n_U\), and \(T\). Those six IVODEs can then be written as derivative expressions either by algebraic manipulation or using a mass matrix.

  3. Finally, in this assignment the yield was optimized with respect to an IVODE variable (the reaction time). When that is the case the reactor design equations only need to be solved once, so long as the interval the variables span includes the optimum value of the variable being optimized. When a reactor response is optimized with respect to something other than an independent or dependent ODE variable, the design equations typically need to be solved multiple times. For example, if the yield in this assignment had been optimized with respect to the heat transfer area, the design equations would need to be solved for many different values of the heat transfer area in order to find one that maximizes the yield.

9.6.4 Optimization of Net Production Rate

The rate expression for liquid-phase reaction (1) is given in equation (2). The rate coefficient displays Arrhenius temperature dependence with a pre-exponential factor of 2.59 x 109 min-1 and an activation energy of 16.5 kcal mol-1. The heat of reaction (1) may be taken to be constant and equal to -22,200 cal mol-1. A solution containing only A at a concentration of 2 M and a temperature of 23 °C is going to be processed in a BSTR. The heat capacity of the solution is approximately constant and equal to 440 cal L-1 K-1, and its density is constant.

\[ A \rightarrow Z \tag{1} \]

\[ r_1 = k_1C_A \tag{2} \]

The 4.0 L BSTR has a perfectly mixed jacket with a volume of 0.5 L, a heat transfer area of 0.6 ft2, and a heat transfer coefficient of 1.13 x 104 cal ft-2 h-1 K-1. Cooling water at 20 °C can be fed to the jacket. The water may be taken to have a constant density of 1 g cm-3 and a constant heat capacity of 1 cal g-1 K-1. A heating coil with a heat transfer coefficient of 3.8 x 104 cal ft-2 h-1 K-1 and a heat transfer area of 0.23 ft2 can be submerged in and extracted from the reacting solution. Saturated steam at 120 °C can admitted to the coil.

The BSTR is charged with 4 L of a 2 M solution of A at 23 °C. Initially there is no flow to the jacket, but it is filled with water, also at 23 °C. To start the batch process, the steam is admitted to the coil, and the coil is immersed in the solution. When the reacting fluid reaches 50 °C, the coil is extracted from the reacting fluid and cooling water flow to the jacket is started. When the reacting fluid reaches a temperature of 25 °C the cooling water flow is stopped, the reactor is drained and preparations for processing the next batch begin. The BSTR pressure is constant throughout all stages of processing. If the turnaround time for the reactor is 25 min, what coolant flow rate will maximize the net rate of production of Z? Plot the conversion of A and the reacting fluid temperature vs. reaction time corresponding to the coolant flow rate that maximizes the net rate and comment on the results.


This is an isolated reactor analysis problem: the reactor is the only equipment involved, the reactions and rate expressions are known, and the questions are about reactor inputs and outputs. To summarize the assignment, I’ll assign appropriate variable symbols to every quantity mentioned in the assignment. This is a batch reactor with 2-stage operating protocol. I’ll use a subscripted “0” to denote initial conditions, a subscripted “1” to denote values at the end of the first stage in the protocol and a subscripted “f” to denote final values. Since there is only one reaction, I won’t use subscripts to index reactions.

This reactor has both a coil for heating and a jacket for cooling. I’ll use a subscripted “ex” to denote properties associated with the cooling jacket and a sub-scripted “coil” to denote properties associated with the heating coil.

9.6.4.1 Assignment Summary

Reaction:

\[ A \rightarrow Z \tag{1} \]

Rate Expression:

\[ r = kC_A \tag{2} \]

Reactor System: BSTR with stage 1 heating and stage 2 cooling.

Reactor Schematic:

Figure 9.10: Schematic representation of a closed (batch), stirred tank reactor with a heating coil (yellow) and a cooling jacket (green).

Quantities of Interest: \(\dot{m}_{ex,opt} = \underset{\dot{m}_{ex}}{\arg\max} \left(r_{Z,net}\right)\), \(\underline{f}_A\big\vert_{\dot{m}_{ex,opt}}\) vs. \(\underline{t}\), and \(\underline{T}\big\vert_{\dot{m}_{ex,opt}}\) vs. \(\underline{t}\) as graphs.

Given and Known Constants: \(k_0\) = 2.59 x 109 min-1, \(E\) = 16.5 kcal mol-1, \(\Delta H\) = -22,200 cal mol-1, \(C_{A,0}\) = 2 M, \(T_0\) = 23 °C, \(\breve{C}_p\) = 440 cal L-1 K-1, \(V\) = 4.0 L, \(V_{ex}\) = 0.5 L, \(A_{ex}\) = 0.6 ft2, \(U_{ex}\) = 1.13 x 104 cal ft-2 h-1 K-1, \(T_{ex,in}\) = 20 °C, \(\rho_{ex}\) = 1 g cm-3, \(\tilde{C}_{p,ex}\) = 1 cal g-1 K-1, \(U_{coil}\) = 3.8 x 104 cal ft-2 h-1 K-1, \(A_{coil}\) = 0.23 ft2, \(T_{coil}\) = 120 °C, \(T_{ex,out,0}\) = 23 °C, \(T_1\) = 50 °C, \(T_f\) = 25 °C, \(t_{turn}\) = 25 min, and \(R\) = 1.987 cal mol-1 K-1.

9.6.4.2 Generation of the Model Equations

The first thing I need to do is to generate the reactor design equations. The operational protocol for this reactor has two stages, so I need to write and simplify the reactor design equations needed to model each stage. The reactor design equations always include at least one mole balance. The general BSTR mole balance is given in Equation 6.8. I’ll write a mole balance for both of the reagents in this system, noting that the sum reduces to a single term since only one reaction is taking place.

\[ \frac{dn_i}{dt} = V \cancelto{\nu_{i,1}r_1}{\sum_j \nu_{i,j}r_j} \]

The BSTR is not isothermal, so I must include an energy balance on the reacting fluid among the design equations. The general BSTR mole balance is given in Equation 6.9. Assuming the liquid to be incompressible and the reactor walls to be rigid, both the pressure and volume will be constant, so their time-derivatives will equal zero. The work associated with mixing the reactor can also be assumed to be negligible. There is only one reaction occurring, so the final sum will reduce to a single term. In addition, the assignment provides the volumetric heat capacity of the entire solution, so the sensible heat term can be written in terms of that heat capacity, Equation 6.10. During the first stage of the operating protocol, heat is exchanged with both the cooling water and the steam, so the rate of heat exchange must be split into two terms. During the second stage the coil is not present.

\[ \cancelto{V \breve{C}_p}{\left(\sum_i n_i \hat C_{p,i} \right)} \frac{dT}{dt} - V\cancelto{0}{\frac{dP}{dt}} - P \cancelto{0}{\frac{dV}{dt}} = \cancelto{\dot{Q}_{ex} + \dot{Q}_{coil}}{\dot Q} - \cancelto{0}{\dot W} - V \cancelto{r_1 \Delta H_1}{\sum_j \left(r_j \Delta H_j \right)} \]

\[ V \breve{C}_p \frac{dT}{dt} = \dot{Q}_{ex} + \dot{Q}_{coil} - V r_1 \Delta H_1 \]

There are two heat exchange fluids in this system. The cooling water in the jacket exchanges sensible heat so the energy balance on it is given by Equation 6.2. During the first stage of the operating protocol, the coolant is not flowing, so \(\dot{m}_{ex} = 0\). During the second stage, the cooling water is flowing. Noting that the heat capacity is constant, the integral is easily evaluated for the second stage.

\[ \rho_{ex} V_{ex} \tilde C_{p,ex}\frac{dT_{ex}}{dt} = -\dot{Q}_{ex} - \dot m_{ex} \int_{T_{ex,in}}^{T_{ex}} \tilde C_{p,ex}dT \]

\[ \rho_{ex} V_{ex} \tilde C_{p,ex}\frac{dT_{ex,out}}{dt} = -\dot{Q}_{ex} - \dot m_{ex}\tilde C_{p,ex} \left(T_{ex,out}-T_{ex,in}\right) \]

The steam in the coil exchanges latent heat, so its temperature is known and constant. As a consequence, the mole balances, energy balance on the reacting fluid and energy balance on the cooling water can be solved independently of the energy balance on the steam. The problem does not ask any questions about the steam flow rate or how much of it condenses, so an energy balance on the steam is not needed.

Momentum balances are not used for stirred tank reactors, so the full set of reactor design equations consists of mole balances on A and Z, an energy balance on the reacting fluid and an energy balance on the cooling water. The mole balances are the same for both stages of operation, but both energy balances change from one stage to the other. For each stage there are four IVODEs and four dependent variables, so it is not necessary to add an IVODE or eliminate a dependent variable.

BSTR Design Equations

Mole balances on A and Z, equations (3) and (4) are the same in heating stage of the protocol and in the cooling stage. An energy balance on the reacting fluid and an energy balance on the coolant in the shell during the heating stage are shown in equations (5) and (6). During the cooling stage of the operating protocol, an energy balance on the reacting fluid and an energy balance on the coolant take the forms shown in equations (7) and (8).

\[ \frac{dn_A}{dt} = -Vr \tag{3} \]

\[ \frac{dn_Z}{dt} = Vr \tag{4} \]

\[ \frac{dT}{dt} = \frac{\dot{Q}_{ex} + \dot{Q}_{coil} - Vr \Delta H}{V \breve C_p} \tag{5} \]

\[ \frac{dT_{ex}}{dt} = -\left(\frac{\dot{Q}_{ex}}{\rho_{ex} V_{ex} \tilde C_{p,ex}}\right) \tag{6} \]

\[ \frac{dT}{dt} = \frac{\dot{Q}_{ex} - Vr \Delta H}{V \breve C_p} \tag{7} \]

\[ \frac{dT_{ex,out}}{dt} = -\left(\frac{\dot{Q}_{ex} + \dot m_{ex} \tilde{C}_{p,ex} \left(T_{ex,out} - T_{ex,in}\right)}{\rho_{ex} V_{ex} \tilde C_{p,ex}}\right) \tag{8} \]

The BSTR design equations are ODEs, so initial values and a stopping criterion are needed to solve them. For this system, two sets of initial values and stopping criteria will be needed because the ODEs need to be solved for each of the two operational stages.

The instant the solution is added to the reactor can be defined as \(t=0\). This marks the start of the heating phase. The molar amount of A, \(n_{A,0}\), can be calculated directly using its known initial concentration and the fluid volume. Only A is present initially, so the initial amount of Z is zero. The reacting fluid temperature, \(T_0\), and the cooling water temperature, \(T_{ex,0}\), at that time are also known. The heating stage ends when the reacting fluid reaches \(T_1\), which is also known. All of the initial values for the first stage are known or can be calculated directly, and only one final value is known, so the BSTR design equations for the first stage are IVODEs.

The cooling phase begins the instant the heating phase ends, so the initial values for the cooling phase are the final values from the heating phase. Assuming that the BSTR design equations for the first stage will be solved before those for the second stage, the initial values for the second stage are all known. The cooling phase ends when the reacting fluid cools down to \(T_f\), so that is the stopping criterion for the cooling phase. Again, the BSTR design equations for the second stage are IVODEs.

Initial Values and Stopping Criteria

Table 9.4: Initial values and stopping criterion for solving the design equations, (3) through (6), for the heating stage of the operating protocol.
Variable Initial Value Stopping Criterion
\(t\) \(0\)
\(n_A\) \(n_{A,0}\)
\(n_Z\) 0
\(T\) \(T_0\) \(T_1\)
\(T_{ex}\) \(T_{ex,0}\)

 

Table 9.5: Initial values and stopping criterion for solving the design equations, (3), (4), (7) and (8), for the cooling stage of the operating protocol.
Variable Initial Value Stopping Criterion
\(t\) \(t_1\)
\(n_A\) \(n_{A,1}\)
\(n_Z\) \(n_{Z,1}\)
\(T\) \(T_1\) \(T_f\)
\(T_{ex,out}\) \(T_{ex,1}\)

Supporting Equations for Calculating the Initial Values

\[ n_{A,0} = C_{A,0}V \tag{9} \]

In order to solve the IVODEs numerically I’ll need write a derivatives function for each of the two operational stages. Each derivatives function will receive \(t\), \(n_A\), \(n_Z\), \(T\), and \(T_{ex,out}\) at the start of an integrations step. For the heating stage, it must use equations (3) through (6) to evaluate the derivatives, and then return them. Before the derivatives can be evaluated, any additional unknowns appearing in those equations must be calculated. Examination of equations (3) through (6) shows that for the heating stage, \(r\), \(\dot{Q}_{ex}\), and \(\dot{Q}_{coil}\) are additional unknowns that must be calculated.

The heat exchange rates can be calculated using the known heat transfer coefficients and heat transfer areas, Equation E.4. The reaction rate can be calculated using the given rate expression, but doing so introduces two more unknowns, \(k\) and \(C_A\). The rate coefficient can be calculated using the Arrhenius expression, Equation 4.8, since the pre-expoenetial factor and activation energy are known. The concentration of reagent A can be calculated using the defining equation for concentration in a closed system, Equation 2.7.

For the cooling stage equations (3), (4), (7), and (8) will be used to evaluate the derivatives. Examination of those equations shows that in addition to \(r\) and \(\dot{Q}_{coil}\), the coolant flow rate, \(\dot{m}_{ex}\), also is an unknown. The reaction rate and the rate of heat transfer to the coil can be calculated as before, but \(\dot{m}_{ex}\) cannot be calculated directly. As discussed in Appendix F.5.2, it will need to be made available to the derivatives function by some means other than passing it as an argument.

Ancillary Equations for Evaluating the Derivatives

During both stages of the operational protocol, the reaction rate can be calculated using the rate expression, equation (2). The rate coefficient and concentration of A appearing in the rate expression can be calculated using equations (10) and (11). The rate of heat transfer to the jacket for both stages can be calculated using equation (12). (Note that Equation (12) uses \(T_{ex,out}\) for the temperature in the jacket, but during the heating stage the coolant is not flowing, and the “out” can be omitted as was done in equation (6) and Tables 9.4 and 9.5.) The coil is only used during the heating stage, and its rate of heat transfer can be calculated using equation (13).

\[ k_1 = k_{0,1}\exp{\left( \frac{-E_1}{RT} \right)} \tag{10} \]

\[ C_A = \frac{n_A}{V} \tag{11} \]

\[ \dot{Q}_{ex} = U_{ex} A_{ex} \left(T_{ex,out} - T\right) \tag{12} \]

\[ \dot{Q}_{coil} = U_{coil} A_{coil} \left(T_{coil} - T\right) \tag{13} \]

At this point, the design equations for the first stage can be solved to find corresponding sets of values of \(t\), \(n_A\), \(n_Z\), \(T\), and \(T_{ex}\) that span the range from their first-stage initial values to their first-stage final values. Then, given a value of \(\dot{m}_{ex}\), the design equations for the second stage can be solved to find corresponding sets of values of \(t\), \(n_A\), \(n_Z\), \(T\), and \(T_{ex}\) that span the range from their second-stage initial values to their second-stage final values. Using the final molar amount of Z, the net rate of production of Z can then be calculated from its definition, Equation 4.2. Similarly, using the molar amount of A, the conversion can be calculated at each time using its definition, Equation 3.4. If the calculations are repeated using a range of values of \(\dot{m}_{ex}\) that is sufficiently wide that \(r_{Z,net}\) passes through a maximum, the maximum net rate and the optimum coolant flow rate can be identified.

9.6.4.3 Calculating the Quantities of Interest

The coolant flow rate that maximizes the net rate of production of Z is one of the quantities of interest. To identify that optimum coolant flow rate, a set of 100 values of \(\dot{m}_{ex}\) between 100 and 250 g min-1 can be created, equations (14) and (15).

\[ N=100 \tag{14} \]

\[ \underline{\dot{m}}_{ex} = \left[\dot{m}_{ex,1},\cdots,\dot{m}_{ex,N}\right]:\,\dot{m}_{ex,1}=100, \dot{m}_{ex,i} \lt \dot{m}_{ex,i+1} \text{ and } \dot{m}_{ex,N} = 250 \tag{15} \]

For each coolant flow rate the design equations for the heating stage can be solved numerically, and the final values of the independent and dependent variables can the be used to define the initial values for the cooling stage. The design equations for the cooling stage can then be solved to get corresponding sets of values of \(t\), \(n_A\), \(n_Z\), \(T\), and \(T_{ex,out}\) during the second, cooling stage. Using the final values, corresponding values of the net rate and the conversion can be calculated using equation (16).

\[ \forall \, \dot{m}_{ex,i} \quad \begin{pmatrix}\displaystyle r_{Z,net,i} = \frac{n_{Z,f,i}}{t_{f,i} + t_{turn}}\\ \\ \displaystyle f_{A,i} = \frac{n_{A,0} - n_{A,f,i}}{n_{A,0}} \end{pmatrix} \tag{16} \]

The resulting set of values for the net rate can be used to identify the coolant flow rate that maximizes the net rate, along with the maximum net rate, equations (17), (18), and (19).

\[ i_{opt} = \underset{i}{\arg\max} \left(r_{Z,net,i}\right) \tag{17} \]

\[ \dot{m}_{ex,opt} = \dot{m}_{ex,i_{opt}} \tag{18} \]

\[ r_{Z,net,max} = r_{Z,net,i_{opt}} \tag{19} \]

If the net rate is maximized at a coolant flow rate of 100 g min-1, the calculations must be repeated using a range of smaller flow rates, and if it is maximized at a coolant flow rate of 250 g min-1 they must be repeated using a range of larger flow rates. Otherwise, the BSTR design equations for both stages can be solved using the optimum coolant flow rate, and the results can be combined to generate corresponding sets of values of \(t\), \(n_A\), \(n_Z\), \(T\), and \(T_{ex,out}\) that span both stages of operation. The corresponding conversion can be calculated using equation (20), and plots of \(\underline{f}_A\big\vert_{\dot{m}_{ex,opt}}\) vs. \(\underline{t}\), and \(\underline{T}\big\vert_{\dot{m}_{ex,opt}}\) vs. \(\underline{t}\) can be generated.

\[ f_{A,i} = \frac{n_{A,0} - n_{A,f,i}}{n_{A,0}} \tag{20} \]

9.6.4.4 Formulation of the Calculations

The BSTR design equations for this assignment are IVODEs, and they can be solved numerically as described in Appendix F.5.2. To complete the assignment as described in the Calculating the Quantities of Interest section, five computer functions can be written.

  • The derivatives function for the heating stage must be passed the independent and dependent variables at the start of an integration step, and it must return the derivatives calculated using equations (3) through (6). To do so, the algorithm must calculate the rate, \(r\), and the two heat transfer rates, \(\dot{Q}_{ex}\) and \(\dot{Q}_{coil}\) and then evaluate the derivatives.

  • The derivatives function for the cooling stage must also be must be passed the independent and dependent variables at the start of an integration step and must return the derivatives, but it must calculate them using equations (3), (4), (7), and (8). To do so, the coolant flow rate must be available and the algorithm must calculate the rate, \(r\), and the heat transfer rate, \(\dot{Q}_{ex}\), before evaluating the derivatives.

  • The BSTR function for the heating stage does not need to be passed any arguments, and it returns corresponding sets of values of the independent and dependent variables spanning the duration of that stage. To do so the algorithm must define the initial values and stopping criterion for the first stage and then call an IVODE solver.

  • The BSTR function for the cooling stage must be passed the final values of the independent and dependent variables in the first (heating) stage and the coolant flow rate as arguments, and it returns corresponding sets of values of the independent and dependent variables spanning the duration of the cooling stage. To do so the algorithm must make the coolant flow rate available to the second stage derivatives function, define the initial values and stopping criterion, and call an IVODE solver.

  • The deliverables function will take no arguments, and it will not return anything. It must then calculate the quantities of interest as described above.

Specifications can be written for each of these functions in preparation for writing them. The function specifications given here assume that the known constants are available wherever they are needed.

Derivatives Function for the Heating Stage Specifications

Arguments: \(t\), \(n_A\), \(n_Z\), \(T\), and \(T_{ex,out}\)

Returns: \(\frac{dn_A}{dt}\), \(\frac{dn_Z}{dt}\), \(\frac{dT}{dt}\), and \(\frac{dT_{ex,out}}{dt}\)

Algorithm:

\[ k_1 = k_{0,1}\exp{\left( \frac{-E_1}{RT} \right)} \tag{10} \]

\[ C_A = \frac{n_A}{V} \tag{11} \]

\[ \dot{Q}_{ex} = U_{ex} A_{ex} \left(T_{ex,out} - T\right) \tag{12} \]

\[ \dot{Q}_{coil} = U_{coil} A_{coil} \left(T_{coil} - T\right) \tag{13} \]

\[ \frac{dn_A}{dt} = -Vr \tag{3} \]

\[ \frac{dn_Z}{dt} = Vr \tag{4} \]

\[ \frac{dT}{dt} = \frac{\dot{Q}_{ex} + \dot{Q}_{coil} - Vr \Delta H}{V \breve C_p} \tag{5} \]

\[ \frac{dT_{ex}}{dt} = -\left(\frac{\dot{Q}_{ex}}{\rho_{ex} V_{ex} \tilde C_{p,ex}}\right) \tag{6} \]

Derivatives Function for the Cooling Stage Specifications

Arguments: \(t\), \(n_A\), \(n_Z\), \(T\), and \(T_{ex,out}\)

Must be Available: \(\dot{m}_{ex}\)

Returns: \(\frac{dn_A}{dt}\), \(\frac{dn_Z}{dt}\), \(\frac{dT}{dt}\), and \(\frac{dT_{ex,out}}{dt}\)

Algorithm:

\[ k_1 = k_{0,1}\exp{\left( \frac{-E_1}{RT} \right)} \tag{10} \]

\[ C_A = \frac{n_A}{V} \tag{11} \]

\[ \dot{Q}_{ex} = U_{ex} A_{ex} \left(T_{ex,out} - T\right) \tag{12} \]

\[ \frac{dn_A}{dt} = -Vr \tag{3} \]

\[ \frac{dn_Z}{dt} = Vr \tag{4} \]

\[ \frac{dT}{dt} = \frac{\dot{Q}_{ex} - Vr \Delta H}{V \breve C_p} \tag{7} \]

\[ \frac{dT_{ex,out}}{dt} = -\left(\frac{\dot{Q}_{ex} + \dot m_{ex} \tilde{C}_{p,ex} \left(T_{ex,out} - T_{ex,in}\right)}{\rho_{ex} V_{ex} \tilde C_{p,ex}}\right) \tag{8} \]

BSTR Function for the Heating Stage Specifications

Arguments: none.

Returns: \(\underline{t}\), \(\underline{n}_A\), \(\underline{n}_Z\), \(\underline{T}\), and \(\underline{T}_{ex}\)

Algorithm:

\[ t_{initial} = 0 \tag{20} \]

\[ n_{A,initial} = n_{A,0} = C_{A,0}V \tag{9} \]

\[ n_{Z,initial} = 0,\, T_{initial} = T_0,\, T_{ex,initial} = T_{ex,0},\, T_{final} = T_1\tag{21} \]

 

\[ \begin{matrix}t_{initial}, n_{A,initial}, n_{Z,initial},\\ T_{initial}, T_{ex,initial}, T_{final},\\ \text{ derivatives function for the heating stage}\\ \Downarrow\\ \text{IVODE solver}\\ \Downarrow\\ \underline{t}, \underline{n}_A, \underline{n}_Z, \underline{T}, \underline{T}_{ex}\\ \end{matrix} \tag{22} \]

Derivatives Function for the Cooling Stage Specifications

Arguments: \(\dot{m}_{ex}\), \(t_{initial}\), \(n_{A,initial}\), \(n_{Z,initial}\), \(T_{initial}\), and \(T_{ex,out,initial}\)

Returns: \(\underline{t}\), \(\underline{n}_A\), \(\underline{n}_Z\), \(\underline{T}\), and \(\underline{T}_{ex,out}\)

Algorithm:

\[ \dot{m}_{ex} \, \Rightarrow \, \text{available to the derivatives function for the cooling stage}\tag{23} \]

 

\[ \begin{matrix}t_{initial}, n_{A,initial}, n_{Z,initial},\\ T_{initial}, T_{ex,out,initial}, T_{final},\\ \text{ derivatives function for the cooling stage}\\ \Downarrow\\ \text{IVODE solver}\\ \Downarrow\\ \underline{t}, \underline{n}_A, \underline{n}_Z, \underline{T}, \underline{T}_{ex}\\ \end{matrix} \tag{24} \]

Deliverables Function Specifications

Arguments: none

Returns: nothing

Algorithm:

\[ \begin{matrix} \text{BSTR function for the Heating Stage}\\ \Downarrow\\ \underline{t}_{heating}, \underline{n}_{A,heating}, \underline{n}_{Z,heating}, \underline{T}_{heating}, \underline{T}_{ex,heating} \end{matrix} \tag{25} \]

 

\[ N=100 \tag{14} \]

\[ \underline{\dot{m}}_{ex} = \left[\dot{m}_{ex,1},\cdots,\dot{m}_{ex,N}\right]:\,\dot{m}_{ex,1}=100, \dot{m}_{ex,i} \lt \dot{m}_{ex,i+1} \text{ and } \dot{m}_{ex,N} = 250 \tag{15} \]

 

\[ \begin{matrix} t_{initial} = \underline{t}_{heating}\big\vert_{T_1},\, n_{A,initial} = \underline{n}_{A,heating}\big\vert_{T_1},\, n_{Z,initial} = \underline{n}_{Z,heating}\big\vert_{T_1} \\ T_{initial} = \underline{T}_{heating}\big\vert_{T_1},\, T_{ex,out,initial} = \underline{T}_{ex,heating}\big\vert_{T_1},\, T_{final} = T_f \end{matrix}\tag{26} \]

 

\[ \forall \, \dot{m}_{ex,i} \quad \begin{pmatrix} \begin{matrix} \dot{m}_{ex,i},\, t_{initial},\, n_{A,initial},\, n_{Z,initial},\\ T_{initial},\, T_{ex,out,initial},\, T_{final} \\ \Downarrow\\ \text{BSTR function for the Cooling Stage}\\ \Downarrow\\ \underline{t}, \underline{n}_A, \underline{n}_Z, \underline{T}, \underline{T}_{ex} \end{matrix}\\ \\ \displaystyle r_{Z,net,i} = \frac{n_{Z,f,i}}{t_{f,i} + t_{turn}} \end{pmatrix} \tag{16} \]

 

\[ i_{opt} = \underset{i}{\arg\max} \left(r_{Z,net,i}\right) \tag{17} \]

\[ \dot{m}_{ex,opt} = \dot{m}_{ex,i_{opt}} \tag{18} \]

 

\[ \begin{matrix} \dot{m}_{ex,opt},\, r_{Z,net,i_{opt}}\\ \Downarrow\\ \text{functions to display and save} \end{matrix} \tag{27} \]

 

\[ \begin{matrix} \dot{m}_{ex,opt},\, t_{initial},\, n_{A,initial},\, n_{Z,initial},\\ T_{initial},\, T_{ex,out,initial},\, T_{final} \\ \Downarrow\\ \text{BSTR function for the Cooling Stage}\\ \Downarrow\\ \underline{t}_{cooling}, \underline{n}_{A,cooling}, \underline{n}_{Z,cooling}, \underline{T}_{cooling}, \underline{T}_{ex,out,cooling} \end{matrix} \tag{28} \]

 

\[ \underline{t} = \left[ \underline{t}_{heating}, \underline{t}_{cooling}\right] \tag{29} \]

\[ \underline{n}_A = \left[ \underline{n}_{A,heating}, \underline{n}_{A,cooling}\right] \tag{30} \]

\[ \underline{f}_A = \frac{n_{A,0} - \underline{n}_A}{n_{A,0}} \tag{31} \]

\[ \underline{T} = \left[ \underline{T}_{heating}, \underline{T}_{cooling}\right] \tag{32} \]

 

\[ \begin{matrix} \underline{t},\, \underline{f}_A,\, \underline{T}\\ \Downarrow\\ \text{plotting functions}\\ \Downarrow\\ \text{requested graphs} \end{matrix} \tag{33} \]

9.6.4.5 Implementation and Execution of the Calculations

The five functions above were written as specified using an appropriate mathematics software package. The calculations were performed by calling the deliverables function. The optimum coolant flow rate was between 100 and 250 g min-1, so it was not necessary to repeat the calculations using a different range of coolant flow rates.

9.6.4.6 Results and Discussion

The net rate of production of Z is plotted as a function of the cooling water flow rate in Figure 9.11. The maximum net rate, 0.063 mol min-1, occurs at a coolant flow rate of 183 g min-1. The conversion and reacting fluid temperature profiles at that coolant flow rate are shown in Figure 9.12 and Figure 9.13. Knowing the initial concentration of A and the Arrhenius parameters, the conversion and temperature profiles were used to calculate the instantaneous reaction rate profile shown in Figure 9.14. It shows that the rate is close to zero when the reactant is charged to the BSTR.

Figure 9.11: Net rate of production of Z as cooling water flow rate varies.
Figure 9.12: Conversion profile during processing with the optimum coolant flow rate.
Figure 9.13: Temperature profile during processing with the optimum coolant flow rate.
Figure 9.14: Instantaneous reaction rate profile during processing with the optimum coolant flow rate.

Four interesting features of Figures 9.12 through 9.14 can be rationalized quatitatively. First, Figure 9.13 shows a kink at approximately 2.5 min. The temperature at that point is 50 °C, and that is the point where the heating stage ends and the cooling stage begins. Just prior to the kink, the temperature was rising rapidly and nearly linearly. That indicates that the temperature rise was primarily due to the heat from the coil. After the kink, the temperature continues to rise, but not as rapidly, because the rate has become large enough that heat is being released by reaction faster than it is being removed by the coolant. The kink in Figure 9.14 occurs for the same reason.

The second interesting observation is that the instantaneous rate (Figure 9.14) reaches a maximum before the temperature (Figure 9.13) does. The reason for the maximum in the instantaneous rate is similar to Example 9.6.1. Before the maximum rate is attained, two opposing factors are affecting the change in the reaction rate. The reactant concentrations are decreasing, which tends to decrease the rate, while the temperature is increasing, which tends to increase the rate. Before the maximum is reached, the temperature effect is greater than the concentration effect. When the two effects become equal, the maximum rate occurs.

In Example 9.6.1 the temperature rose continually because the reactor was adiabatic, but here Figure 9.13 shows that the temperature does eventually reach a maximum and then decrease. This is due to heat removal by the cooling jacket. The third interesting feature of the graphs is that the temperature continues to rise for a short period of time after the rate has reached its maximum value. Even though the rate has begun to decrease due to depletion of the reactants, the reaction is still releasing heat faster than the cooling fluid can remove it. When the temperature maximum is reached, the release of heat by reaction has become equal to the removal of heat by the coolant. After the maximum, the coolant removes heat faster than it is generated, so that after ca. 100 min the temperature reaches 25 °C and processing ends.

The conversion at this point is over 99%, but the fourth interesting feature is that the conversion is well over 90% after approximately 30 min. Put differently approximately 90% conversion occurs in the first 30 minutes, but it only increases 9% in the final 70 minutes. This suggests that it might be possible to optimize the system further.

In this assignment, the net rate was maximized with respect to the coolant flow rate. The net rate will also be affected by the duration of the heating stage, or equivalently the temperature at which the coil is removed from the reacting fluid and the coolant flow starts. If the temperatures of the cooling water or the steam could be varied, they, too, would be expected to affect the net rate. In a more realistic optimization, it would not be the net rate that would be optimized. Instead, the rate of financial profit would be maximized. The rate of financial profit would depend on the net rate of the reaction, but it would also depend on the selling price of the product, the costs of the reactant, steam, and cooling water, and a number of other factors. That type of multidimensional optimization is beyond the scope of this book, but at its core, it is similar to the simple optimization illustrated in this assignment.

9.7 Symbols Used in Chapter 9

Symbol Meaning
\(i\) index denoting a reagent.
\(\dot{m}_{ex}\) mass flow rate of the heat exchange fluid.
\(n_i\) molar amount of reagent \(i\), an additional subscripted 0 denotes the initial molar amount.
\(r_{i,net}\) net rate of generation of reagent \(i\).
\(r_j\) rate of reaction \(j\) per volume of reacting fluid.
\(t\) time.
\(t_{rxn}\) length of time during which reaction occurs.
\(t_{turn}\) turnaround time.
\(A\) heat transfer area.
\(C_i\) concentration of reagent \(i\).
\(\hat{C}_{p,i}\) molar, constant-pressure heat capacity of reagent \(i\)
\(\tilde{C}_{p,ex}\) mass-specific, constant-pressure heat capacity of the heat exchange fluid.
\(\breve{C}_p\) volume-specific, constant-pressure heat capacity of the reacting fluid, additional subscripts denote fluids other than the reacting fluid.
\(M_{ex}\) molecular weight of the heat exchange fluid.
\(P\) pressure of the reacting fluid, a subscripted \(i\) denotes the partial pressure of reagent \(i\).
\(\dot{Q}\) rate of heat transfer to the reacting fluid, additional subscripts distinguish between different sources of the heat.
\(R\) ideal gas constant.
\(T\) temperature, no subscript denotes the reacting fluid, a subscripted \(ex\) denotes the heat exchange fluid, an additional subscripted \(0\) denotes an initial temperature, \(in\) an inlet temperature and \(f\) a final temperature.
\(U\) heat transfer coefficient.
\(V\) reacting fluid volume, a subscripted \(ex\) denotes the heat exchange fluid.
\(\dot{W}\) rate at which the reacting fluid performs work on its surroundings.
\(\gamma\) fraction of the heat exchange fluid that undergoes phase change.
\(\nu_{i,j}\) stoichiometric coefficient of reagent \(i\) in reaction \(j\).
\(\rho\) density, no subscript denotes the reacting fluid, a subscripted \(ex\) denotes the heat exchange fluid.
\(\Delta H_j\) heat of reaction \(j\).