23  A Segregated Flow Reactor Model

Chapter 22 introduced the cumulative and differential age distribution functions. It showed how to measure the cumulative age distribution function, \(F\), and how it is related to the differential age distribution function, \(dF\). It also mentioned that the age distribution functions can be used to construct models for non-ideal reactors. This chapter describes one way to model a non-ideal flow reactor using its cumulative age distribution function. It is one variant of a group of segregated flow reactor (SFR) models. Only steady-state non-ideal reactors that operate at constant pressure and either adiabatically or isothermally are considered. This is not as limiting as it might first seem. Typically stirred tanks operate at constant pressure, and some tubular reactors operate with negligible pressure drop.

23.1 SFR Model Assumptions

The SFR model makes use of the concept of fluid elements that was first introduced in Chapter 22. The fluid entering the reactor is assumed to be segregated into many, many small volumes called fluid elements. Fluid elements aren’t real, observable things; they are hypothetical constructs. As such, the volume of a fluid element is not known, beyond that it is very small. However, because fluid elements are formed by splitting the non-ideal reactor feed into small parcels, the initial concentrations of the reagents and the initial temperature are known. They are equal to the concentrations and temperature of the feed to the non-ideal reactor.

The SFR model assumes that while the fluid elements are within the reactor, they remain perfectly mixed internally, but they do not mix with any other fluid elements. The mixing within fluid elements is sometimes called micro-mixing while mixing between different fluid elemnts is macro-mixing.

The SFR model makes three additional assumptions. The first is that the fluid elements have different residence times. They do not all spend the same amount of time within the reactor. The second additional assumption is that the residence times of the fluid elements are assigned according to the age distribution function for the non-ideal reactor. The final assumption is that the fluid elements become perfectly macro-mixed when they leave the reactor.

23.2 Fluid Elements as BSTRs

Each fluid element is a very small BSTR. It is perfectly mixed and no mass enters or leaves because it doesn’t mix with other fluid elements. By choosing the volume of a fluid element as a basis, the ideal BSTR reactor design equations can be solved to determine the composition, temperature, and pressure within a fluid element as a function of its residence time (age). As previously noted, in this chapter the pressure in the reactor is assumed to be constant.

If the non-ideal reactor being modeled operates isothermally and the reacting fluid is an incompressible, ideal liquid mixture, the only BSTR design equations that are needed to model a fluid element are the BSTR mole balances. If the reacting fluid is an ideal gas, a differential form of the ideal gas law may need to be added to the design equations to account for expansion or contraction. If the non-ideal reactor operates adiabatically, an energy balance on the reacting fluid will also be needed.

The analysis of a fluid element using the BSTR design equations is identical to the analysis of actual BSTRs as described in Chapter 9. The design equations will be IVODEs. Having chosen the fluid element volume as a basis, the initial molar amounts of the reagents can be calculated since the initial concentrations and temperature are known. Depending on the whether the reacting fluid is a liquid or a gas and whether the reactor is isothermal or adiabatic, the molar amounts of the reagents, the temperature, and the fluid volume may all change as a function of reaction time. The molar amounts and volume that are found by solving the BSTR design equations are extensive quantities, so they only apply for the chosen basis. However, if they are converted to intensive quantities, they can be used irrespective of the actual size of a fluid element. Hence the molar amounts and volume found by solving the BSTR design equations are typically used to calculate the concentrations of the reagents. Simlarly, if the reacting fluid volume changes as a function of reaction time or age, \(t_r\), it is converted to a relative volume, \(\gamma\), as shown in Equation 23.1.

\[ \gamma = \frac{V\big\vert_{t_r}}{V_{basis}} \tag{23.1}\]

To summarize, even though the volume of a fluid element is not known, the reagent concentrations in it, its temperature, and its relative volume can all be calculated as a function of reaction time, or equivalently, fluid age, by solving the ideal BSTR design equations.

23.3 The Segregated Flow Reactor Model

If \(\dot{V}_{in}\) is the volumetric flow rate into the non-ideal reactor, then \(\dot{V}_{in} dF\big\vert_{\lambda}\) is the inlet volumetric flow rate of fluid that will remain in the reactor for a reaction time between \(\lambda\) and \(\lambda + d\lambda\). The corresponding outlet volumetric flow rate of fluid that was in the reactor for a reaction time between \(\lambda\) and \(\lambda + d\lambda\) is then \(\dot{V}_{feed} \gamma \big\vert_{\lambda} dF\big\vert_{\lambda}\). The total outlet volumetric flow rate is then the sum of the flow rate for all possible reaction times or fluid ages. Since the fluid age intervals are differential, the sum takes the form of an integral, leading to Equation 23.2 for the total volumetric flow rate at the outlet from the non-ideal reactor. (The second form of the equations results from substituting Equation 22.1 into the first form.) Note that for an incompressible liquid, \(\gamma\) would be constant and equal to 1, and because \(\int_0^\infty dF = 1\) (see Chapter 22), \(\dot{V}_{out} = \dot{V}_{in}\), as expected.

\[ \dot{V}_{out} = \dot{V}_{in} \int_0^\infty \gamma dF = \dot{V}_{in} \int_0^\infty \gamma \frac{dF}{d\lambda} d\lambda \tag{23.2}\]

In a similar manner, the outlet molar flow rate of reagent \(i\) can be found by multiplying the outlet volumetric flow rate by the concentration of \(i\) as shown in Equation 23.3.

\[ \dot{n}_{i,out} = \dot{V}_{in} \int_0^\infty \gamma C_i dF = \dot{V}_{in} \int_0^\infty \gamma C_i \frac{dF}{d\lambda} d\lambda \tag{23.3}\]

If the non-ideal reactor operates isothermally, the outlet temperature is equal to the inlet temperature. However if the non-ideal reactor is adiabatic, the outlet temperature must be calculated. From above, the outlet molar flow rate of reagent \(i\) with an age between \(\lambda\) and \(\lambda + d\lambda\) is equal to \(\dot{V}_{in} \gamma \big\vert_{\lambda} C_i \big\vert_{\lambda} \frac{dF}{d\lambda}\bigg\vert_{\lambda} d\lambda\), and it will be at a temperature, \(T\big\vert_{\lambda}\). When the fluid elements macro-mix as they leave the reactor, each of these molar flows must gain or lose sensible heat to bring it to the final temperature, \(T_{out}\), of the fluid leaving the non-ideal reactor. If the macro-mixing is adiabatic, the sum of the sensible heats is equal to zero, Equation 23.4.

\[ \sum_i \dot{V}_{in} \int_0^\infty \gamma C_i \frac{dF}{d\lambda} \left(\int_T^{T_{out}} \hat{C}_{p,i}dT \right) d\lambda = 0 \tag{23.4}\]

In terms of molar heat capacities, this serves as implicit equation for \(T_{out}\) shown in Equation 23.5. For liquid systems, the volumetric heat capacity of the fluid is sometimes used, leading to the implicit equation for \(T_{out}\) shown in Equation 23.6.

\[ T_{out} : \, \dot{V}_{in} \int_0^\infty \gamma \frac{dF}{d\lambda} \left( \sum_i C_i \int_T^{T_{out}} \hat{C}_{p,i} dT \right)d\lambda = 0 \tag{23.5}\]

\[ T_f : \, \dot{V}_{in} \int_0^\infty \gamma \frac{dF}{d\lambda} \left(\int_T^{T_{out}} \breve{C}_pdT \right)d\lambda = 0 \tag{23.6}\]

23.4 Tabular Age Distribution Functions and Numerical Integration

The SFR model is based upon the cumulative age distribution function for the non-ideal reactor being modeled. When the cumulative age distribution function is available as an analytical function, evaluating \(F\) and \(\frac{dF}{d \lambda}\) at any value of \(\lambda\) is straightforward. Often, however, the cumulative age distribution function is measured using stimulus-response experiments (see Chapter 22). In that case, the cumulative age distribution function is typically available as a table with columns containing corresponding values of \(\lambda\) and \(F\).

In this chapter this form of the cumulative age distribution function will be referred to as a table. It should be mentioned that it could also be described as a pair of corresponding vectors. In that case, “adding a column” to the tabular form of the cumulative age distribution function is the same as creating an additional vector that corresponds to the original pair of vectors.

When using a tabular cumulative age distribution function that was generated experimentally, it is important to ensure that it spans a sufficiently large range of ages so that \(F\) ranges from zero to one. When using the SFR model, it is also necessary to add a third column to the table that contains the values of \(\frac{dF}{d\lambda}\) corresponding to each of the ages.

23.4.1 The Tabular Form of the Cumulative Age Distribution Function

By definition, \(F\left(0\right) = 0\), so the first row in the tabular form of the cumulative age distribution always should be \(\lambda\) = 0 and \(F\) = 0. Eventually all of the fluid leaves the reactor, so \(F\) should always equal 1 in the last row of the tabular age distribution. In principle, the age of the very last bit of fluid to leave could equal infinity, but since it is determined experimentally, the largest age in the tabular age function, \(\lambda_f\), will be finite.

That is, \(\lambda_f\) will equal the elapsed time when the last measurement in the stimulus-response experiments was made. Two situations are possible. In one situation, the value of \(F\) will have reached a value of 1 (within experimental error) by that time. By definition, \(F\) increases monotonically as \(\lambda\) increases, so once the maximum value of 1 has been reached, \(F\) will hold constant (to within experimental error) at that value for all later experimental measurements.

However, if the stimulus-response experiments were not continued for a sufficient length of time, the final value of \(F\) in the tabular form of the cumulative age distribution function will be less than 1. In order to use the SFR model, the cumulative age distribution function must range from \(F=0\) to \(F=1\). So if the final value of \(F\) is less than one, a row should be added to the end of the table. In the added row, \(F\) should equal 1, but the proper corresponding value of \(\lambda\) is unknown. One way to choose the value, \(\lambda_f\), to enter in that last row is to linearly extrapolate using the preceding two entries in the experimental cumulative age distribution table to find the value of \(\lambda\) where \(F\) equals 1, Equation 23.7.

\[ \lambda_f = \left(1 - F\big\vert_{\lambda_{f-2}}\right)\frac{\lambda_{f-1} - \lambda_{f-2}}{F\big\vert_{\lambda_{f-1}} - F\big\vert_{\lambda_{f-2}}} + \lambda_{f-2} \tag{23.7}\]

To summarize, the cumulative age distribution function will be a table with two columns. The first column should contain values of \(\lambda\), and the value in the first row of that column should be zero. The second column should contain corresponding values of \(F\). The value in the first row of the \(F\) column should be zero, and the value in the last row of the \(F\) column should equal 1 to within experimental uncertainty.

23.4.2 The Tabular Form of the Differential Age Distribution Function

The derivative of the age function, \(\frac{dF}{d\lambda}\), is needed when using the SFR model. A column representing \(\frac{dF}{d\lambda}\) can be added to the tabular form of the cumulative age distribution function quite easily. For example if forward differences are used to approximate the derivative, Equation 23.8 is used to calculate \(\frac{dF}{d\lambda}\) for all but the last row of the table. Since the age function has reached \(F=1\) in the last row and will remain constant at that value, the derivative in the last row is equal to zero, Equation 23.9

\[ \frac{dF}{d\lambda}\bigg\vert_{\lambda_k} = \frac{F_{k+1} - F_k}{\lambda_{k+1} - \lambda_k} \tag{23.8}\]

\[ \frac{dF}{d\lambda}\bigg\vert_{\lambda_f} = 0.0 \tag{23.9}\]

23.4.3 Numerical Integration

Using the SFR model equations, 23.2 through 23.6, requires the evaluation of an integral over all possible fluid element ages. The quantities in the integrands (\(\gamma\), \(dF\), \(\frac{dF}{d\lambda}\), \(C_i\), and \(T\)) vary with fluid age. Even if the age distribution is avialable as an analytic function, numerical integration will likely be required. That is so because the fluid element design equations can only be solved analytically for a few, very simple systems. Hence, analytical expressions for \(\gamma\), \(C_i\), and \(T\) as functions of fluid element age will not be available.

The integrals indicate infinity as the upper limit, but the critical factor is that \(F\) should equal 1.0 at the upper limit of the integral. In other words, the integral must include all fluid elements. As already mentioned, if the age function was generated experimentally there will be a finite age, \(\lambda_f\), where \(F\) becomes equal to 1.0 to within experimental error. If an analytical age function is available, the upper limit should be set to a finite value, \(\lambda_f\), where \(F\) is essentially equal to 1.0. This is necessary because the finite element design equations will be solved numerically using the fluid element age as the stopping criterion. If the upper limit is left as infinity, the numerical solution of the design equations likely will fail. So, irrespective of whether an analytic or tabular cumulative age distribution is available, the upper limit of the integrals in Equations 23.2 through 23.6 should be changed to \(\lambda_f\).

The trapezoid rule can be used to perform the integration numerically. Doing so requires that the integrand be evaluated at a number of points starting at \(\lambda = 0\) and ending at \(\lambda = \lambda_f\). Here these will be called integrand evaluation points. If there are \(N_{\lambda}\) integrand evaluation points, and if \(y_k\) representes the value of the integrand at \(\lambda_k\), the value of the integral is given by Equation 23.10. (It can be noted that \(\lambda_f = \lambda_{N_\lambda}\)).

\[ \int_0^{\lambda_f} y d\lambda = 0.5 \sum_{k=2}^{N_\lambda} \left[ \left(y_k + y_{k-1}\right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{23.10}\]

It remains to choose the integrand evaluation points. If the age distribution is available as an analytical function, the engineer can decide how many integrand evaluation points to use. Generally, the greater the number of integrand evaluation points, the greater the accuracy. The integrand evaluation points do not need to be evenly spaced. If \(N_\lambda\) integrand evaluation points are selected, it is good practice to repeat the analysis using \(2N_\lambda\) integrand evaluation points and check that the result is not significantly affected.

If the age distribution is available in tabular form, it makes sense to use the ages in the table as the integrand evaluation points because \(\frac{dF}{d\lambda}\) is known at those points. If different integrand evaluation points are chose, it will be necessary to interpolate within the table to find \(\frac{dF}{d\lambda}\).

In either case, once the integrand evaluation points have been selected, the finite element design equations will need to be solved to find \(C_i\), \(T\), and \(\gamma\) for fluid element ages equal to each \(\lambda_k\). These values will be needed, along with \(\frac{dF}{d\lambda}\) at each \(\lambda_k\), to evaluate the integrals in the SFR model equations using the trapezoid rule.

23.5 Numerical Implementation of the SFR Model

Computer code to implement the SFR model numerically could be structured in a number of ways. The structure described here recognizes that in the SFR model, the “reactors” are actually fluid elements. As such, a fluid element model is first defined. This is effectively an ideal BSTR model like those introduced in Chapter 9. An SFR model is then defined. It uses the results from the fliud element model to calculate the non-ideal reactor product molar flow rates and temperature.

23.5.1 The Fluid Element Model

The fluid elements are perfectly mixed internally and have no material inputs or outputs (such as from mixing with other fluid elements). That is, they are effectively BSTRs. Here they are assumed to be adiabatic and to operate at constant pressure. As such the fluid element model is virtually the same as those used to analyze actual BSTRs.

The initial concentrations of the reagents in the fluid element and the temperature are equal to those of the feed to the non-ideal flow reactor. Because fluid elements are hypothetical constructs and not real entities, their volume is not known other than it is very small. Thus, the fluid element model requires that a fluid element volume be assumed and used as a basis for the calculations.

Having chosen the volume, the initial molar amounts can be calculated from their concentrations. The BSTR design equations then need to be solved for each integrand evaluation point. That is the equations need to be solved using each \(\lambda_k\) as the stopping criterion for the time.

For any one \(\lambda_k\), solving the BSTR design equations will yield corresponding sets of values of the reaction time, molar amount of each reagent, temperature, and, for gases, volume. The final values, where the time reaches \(\lambda_k\) are of interest here. Because the initial fluid element volume was chosen as a basis, the final molar amounts and volume must be converted to intensive variables before they are returned. Typically they are converted to reagent concentrations and a relative volume.

Once the concentrations, temperature and relative volume have been calculated for every integrand evaluation point, they are returned by the fluid element function. If a tabular age distribution function is being used, the concentrations, temperatures and relative volumes can be added as additional table columns, if desired.

A fluid element function can be written to perfrom these calculations. Since the BSTR design equations are IVODEs and they are solved numerically within the fluid element model, the fluid element model will also include a derivatives function.

23.5.2 The SFR Model

The primary purpose of the SFR model is to use the results from the finite element model to calculate the outlet molar flow rates and temperature for the non-ideal reactor using Equation 23.3 and either Equation 23.5 or Equation 23.6. Examination of those equations shows that before they can be used, the derivative, \(\frac{dF}{d\lambda}\), must be calculated for each integrand evaluation point. This can be accomplished using Equations 23.8 and 23.9.

The integrals in Equations 23.3, 23.5, and 23.6 then can be evaluated using the trapezoid rule. That is, the molar flow rates leaving the non-ideal reactor can be calculated using Equation 23.11 and the outlet temperature can be calculated using either Equation 23.12 or Equation 23.13.

\[ \dot{n}_{i,out} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(\gamma C_i \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( \gamma C_i \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{23.11}\]

\[ T_{out} \,: 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left(\begin{matrix} \left(\gamma \frac{dF}{d\lambda} \left( \sum_i C_i \int_T^{T_{out}} \hat{C}_{p,i} dT \right)\right)\bigg\vert_{\lambda_k} \\+ \left(\gamma \frac{dF}{d\lambda} \left( \sum_i C_i \int_T^{T_{out}} \hat{C}_{p,i} dT \right) \right)\bigg\vert_{\lambda_{k-1}} \end{matrix} \right) \left(\lambda_k - \lambda_{k-1}\right) \right] = 0 \tag{23.12}\]

\[ T_{out} \,: 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left(\begin{matrix} \left( \gamma \frac{dF}{d\lambda} \left(\int_T^{T_{out}} \breve{C}_pdT \right) \right)\bigg\vert_{\lambda_k} \\+ \left( \gamma \frac{dF}{d\lambda} \left(\int_T^{T_{out}} \breve{C}_pdT \right) \right)\bigg\vert_{\lambda_{k-1}} \end{matrix} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] = 0 \tag{23.13}\]

An SFR function can be written to perform these calculations. The equation for the outlet temperature, either Equation 23.12 or Equation 23.13 is implicit. It can be solved for \(T_{out}\) using an ATE solver. Consequently the SFR model also will include a residuals function.

23.6 Advantages and Disadvantages of SFR Models

A significant advantage of the SFR model presented here is that it is based on the cumulative age distribution function for the non-ideal reactor. As seen in this chapter, when the non-ideal reactor operates adiabatically at constant pressure, this results in an SFR model that does not have any adjustable parameters that need to be estimated. The assumptions of adiabatic and constant pressure operation apply to many real, non-ideal reactors. Stirred tanks normally operate at constant pressure, and often for tubular reactors the pressure drop along the length of the reactor is negligible.

A significant disadvantage of the SFR model is that using it to model non-ideal reactors where the pressure drops or heat is exchanged is not straightforward. If the pressure in the steady-state, non-ideal reactor varies with position, an equation for the variation of the pressure with reaction time could be added to the SFR model above. For example, if the non-ideal reactor is tubular and the pressure drops along its length, the pressure vs. length relationship could be introduced to the SFR model as pressure vs. reaction time in the fluid element. The problem with that approach is that the fluid elements have different ages, and that means they will have different final pressures.

Adapting the SFR model to include heat exchange is not as obvious. One option would be to assume that the fluid element in the SFR model has the same heat transfer area to initial fluid element volume as the non-ideal reactor being modeled. There are problems with this approach, too. For example, the fluid element volume may change as the reaction proceeds. In Reaction Engineering Basics the SFR model is only used for isothermal or adiabatic non-ideal reactors that operate at steady-state and constant pressure.

23.7 Examples

The following two examples illustrate the use of SFR models to simulate non-ideal reactors. In Example 23.7.1 the reacting fluid is an incompressible liquid and in Example 23.7.2 it is an ideal gas. In both examples the cumulative age distribution is available in tabular form.

23.7.1 An SFR Model with a Liquid-Phase Reaction

Suppose the non-ideal reactor from Example 22.5.1 is going to be used to convert an aqueous solution of A to Z as indicated in equation (1). The feed rate to the 10 L reactor is 3 L min-1, and the cumulative age distribution found in that example, cum_age_dist_fcn.csv, applies. The feed temperature is 300 K and it contains only reagent A at a concentration of 0.5 M.

The rate expression is given in equation (2) where the pre-exponential factor is equal to 5.29 x 109 L mol-1 min-1 and the activation energy is 12,100 cal mol-1. The heat of reaction (1) is -24,400 cal mol-1. The heat capacity of the reacting fluid may be taken to be constant at 1.0 cal ml-1 K-1, and the density to be constant at 1.0 g ml-1.

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

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

According to the segregated flow reactor model, what conversion of A and outlet temperature can be expected if the reactor operates adiabatically at a pressure of 1 atm? When these results are calculated using the SFR model, the volume of the non-ideal reactor is not used. Explain how this is possible.


The assignment narrative instructs me to use an SFR model. I’ll begin by summarizing the assignment. I’ll use subscripted “feed” and “prod” to indentify values for the non-ideal reactor feed and product streams, and I’ll use an overbar to indicate a set of values for a quantity. Looking ahead, I know that I will need to assume a fluid element volume as a basis, so I will include that in the summary, using a subscripted “fe” to indicate fluid element. While the volume of fluid elements is assumed to be very small, I can use any convenient value as a basis because I’m only going to calculate intensive quantities.

23.7.1.1 Assignment Summary

Reaction:

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

Rate Expression:

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

Reactor System: Steady-state, adiabatic, constant pressure SFR.

Quantities of Interest: \(f_{A,prod}\) and \(T_{prod}\).

Given and Known Constants: \(V_{non-ideal}\) = 10 L, \(\dot{V}_{feed}\) = 3 L min-1, \(T_{feed}\) = 300 K, \(C_{A,feed}\) = 0.5 mol L-1, \(C_{Z,feed}\) = 0.0, \(P\) = 1 atm, \(k_0\) = 5.29 x 109 L mol-1 min-1, \(E\) = 12,100 cal mol-1, \(\Delta H\) = -24,400 cal mol-1, \(\breve{C}_p\) = 1.0 cal ml-1, K-1, \(\rho\) = 1.0 g ml-1, and \(\overline{F}\) vs. \(\overline{\lambda}\).

Basis: \(V_{fe}\) = 1.0 L

23.7.1.2 Mathematical Formulation of the Analysis

The formulation presented here assumes that the given and known constants identified in the assignment summary, above, are available at any point in the analysis.

In the SFR model, the “reactor” is a fluid element. Fluid elements are effectively adiabatic, constant pressure BSTRs, so I’ll begin by generating a fluid element model. The general form of the BSTR mole balance is given in Equation 6.8. In the present system, only one reaction is taking place, so the summation reduces to a single term and it isn’t necessary to index the reaction.

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

The general BSTR energy balance is given in Equation 6.9. The volumetric heat capacity for the reacting fluid is known and can be used in place of the summation over the molar heat capacities. The pressure is constant, so its time-derivative is zero. Assuming the reacting liquid to be incompressible, the derivative of the volume with respect to time is also zero. There are no shafts or pistions, so the work done by the fluid element may be considered to be negligible, and since the reactor is adiabatic, the heat term equals zero. As above, there is only one reaction, so the final sum reduces to a single term where it is not necessary to index the reaction.

\[ \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{0}{\dot{Q}} - \cancelto{0}{\dot{W}} - V \cancelto{r\Delta H}{\sum_j \left(r_j \Delta H_j \right)} \]

The design equations are IVODEs, and they will be solved numerically, so I’ll rearrange the energy balance in the form of a derivative expression.

\[ \frac{dT}{dt} = - \frac{r \Delta H}{\breve{C}_{p}} \]

Fluid Element Model

Design Equations

\[ \frac{dn_A}{dt} = -V_{fe} r \tag{3} \]

\[ \frac{dn_Z}{dt} = V_{fe} r \tag{4} \]

\[ \frac{dT}{dt} = - \frac{r \Delta H}{\breve{C}_{p}} \tag{5} \]

The design equations are IVODEs, and I’ll solve them using an IVODE solver. To do so, I will need initial values, a stopping criterion, and a derivatives function. I can define \(t=0\) as the instant the fluid is placed in the fluid element. The initial values are then equal to the molar amounts of A and Z and the temperature at that time. The initial concentrations of A and Z are equal to their concentrations in the feed to the non-ideal reactor, and the temperature is equal to the non-ideal reactor feed temperature. The initial moles of A are then equal to the fluid element volume times its initial concentration. The initial moles of Z are equal to zero because its concentration in the feed is zero.

The fluid ages, \(\lambda_k\), in the tabular cumulative age distribution will be used as the integrand evaluation points. That is, I will need to solve the design equations using each \(\lambda_k\) as the final value of the time.

It is also important to check that the cumulative age distribution function spans the full range from \(F=0\) to \(F=1\). Examination of the table for this assignment shows that it does.

In addition to the initial values and stopping criterion, I’ll also need to write a derivatives function that I can provide to the IVODE solver. The derivatives function will be called by the IVODE solver, and the solver will assume that the only arguments are values of the independent (\(t\)) and dependent (\(n_A\), \(n_Z\), and \(T\)) variables at the start of an integration step. The solver will further assume that the derivatives function will return only the derivatives of the dependent variables with respect to the independent variable, and nothing else.

Before the derivatives can be evaluated, any unknown quantities appearing in them must be calculated. Looking at the design equations, the only unknown I see is \(r\). The rate can be calculated using the given rate expression. That will require using the Arrhenius expression, Equation 4.8, and the definition of concentration in a closed system, Equation 1.7, first.

Initial Values and Stopping Criteria:

Table 23.1: Initial values and stopping criterion for solving the fluid element design equations, (3) through (6).
Variable Initial Value Stopping Criteria
\(t\) \(0\) each \(\lambda_k\) in \(\overline{\lambda}\)
\(n_A\) \(C_{A,feed}V_{fe}\)
\(n_Z\) \(0.0\)
\(T\) \(T_{feed}\)

Derivatives Function:

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

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

Algorithm:

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

\[ k = k_0 \exp{\left(\frac{-E}{RT}\right)} \tag{7} \]

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

\[ \frac{dn_A}{dt} = -V_{fe} r \tag{3} \]

\[ \frac{dn_Z}{dt} = V_{fe} r \tag{4} \]

\[ \frac{dT}{dt} = - \frac{r \Delta H}{\breve{C}_{p}} \tag{5} \]

At this point I can write a fluid element function that solves the design equations. I’m going to solve them using different values of \(\lambda_k\) as the stopping criterion. However, since the value of \(\lambda_k\) are those in the cumulative age distribution table, and they are known constants, I don’t need to pass them to this function as arguments.

Solving the design equations using any one \(\lambda_k\) as the stopping criterion will yield corresponding sets of values of \(t\), \(n_A\), \(n_Z\), and \(T\) spanning the range from \(t=0\) to \(t=\lambda_k\). I want to return the final values where the reaction time is equal to \(\lambda_k\).

The fluid element volume was chosen as a basis, so I don’t want to return the final molar amounts. Instead, I will divide them by the volume and return the molar concentrations.

I am going to write this function so that it only needs to be called once and it will return the concentrations and temperatures for all of the integrand evaluation points at once. To do that, each time I calculate the concentrations and temperature for one of the ages, I’ll add them to vectors. Then, after all of the calculations are complete, I’ll return those vectors, and they will correspond to the ages in the tabular cumulative age distribution function.

Fluid Element Function:

Arguments: none

Returns: \(\overline{C}_A\), \(\overline{C}_Z\), and \(\overline{T}\)

Algorithm:

Repeat (9) through (12) for each \(\lambda_k\) in \(\overline{\lambda}\).

\[ \begin{matrix} \text{initial values, stopping criterion, derivatives function} \\ \Downarrow \\ \text{IVODE solver} \\ \Downarrow \\ \overline{t}, \overline{n}_A, \overline{n}_Z, \overline{T} \end{matrix} \tag{9} \]

 

\[ C_{A,k} \text{ in } \overline{C}_A = \frac{n_A \big\vert_{t=\lambda_k}}{V_{fe}} \tag{10} \]

\[ C_{Z,k} \text{ in } \overline{C}_Z = \frac{n_Z \big\vert_{t=\lambda_k}}{V_{fe}} \tag{11} \]

\[ T_k \text{ in } \overline{T} = T\big\vert_{t=\lambda_k} \tag{12} \]

The next thing I need is the SFR model. According to it, the non-ideal reactor outlet molar flow rates are given by Equation 23.3. I will be using the trapezoid rule to evaluate the integral, in which case the equation for the outlet molar flow rates takes the form given in Equation 23.11

\[ \dot{n}_{i,out} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(\gamma C_i \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( \gamma C_i \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \]

The reacting fluid is an incompressible liquid, so the fluid element volume will not change during reaction, and \(\gamma\) will be constant and equal to 1.0.

\[ \dot{n}_{i,out} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left( C_i \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( C_i \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \]

I know the volumetric heat capacity, so the non-ideal reactor product temperature is given by Equation 23.13.

\[ T_{out} \,: 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left(\begin{matrix} \left( \gamma \frac{dF}{d\lambda} \left(\int_T^{T_{out}} \breve{C}_pdT \right) \right)\bigg\vert_{\lambda_k} \\+ \left( \gamma \frac{dF}{d\lambda} \left(\int_T^{T_{out}} \breve{C}_pdT \right) \right)\bigg\vert_{\lambda_{k-1}} \end{matrix} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] = 0 \]

As above, \(\gamma\) will be constant and equal to 1.0. The heat capacity is a constant, so it can be factored out of the sensible heat integral which can then be evaluated analytically. In addition, this is an implicit equation for \(T_{out}\), so I will write in in the form of a residual expression.

\[ 0 = \epsilon = 0.5 \dot{V}_{feed} \breve{C}_p \sum_{k=2}^{N_\lambda} \left[ \left( \left( \frac{dF}{d\lambda} \left( T_{out} - T \right) \right)\bigg\vert_{\lambda_k} + \left( \frac{dF}{d\lambda} \left(T_{out} - T \right) \right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] = 0 \]

This is an implicit equation that will be solved using an ATE solver. To do so, I’ll need to provide a guess for the solution and a residuals function. Since the reactor is adiabatic and the reaction is exothermic, I’ll guess that the product temperature is 10 K greater than the feed temperature. If the ATE solver fails to converge, I may need to revise that guess.

Segregated Flow Reactor Model

Model Equations:

\[ \dot{n}_{A,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(C_A \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( C_A \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{13} \]

\[ \dot{n}_{Z,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(C_Z \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( C_Z \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{14} \]

\[ 0 = \epsilon = 0.5 \dot{V}_{feed} \breve{C}_p \sum_{k=2}^{N_\lambda} \left[ \left( \begin{matrix} \left( \frac{dF}{d\lambda} \left(T_{Prod} - T \right)\right)\bigg\vert_{\lambda_k} \\+ \left( \frac{dF}{d\lambda} \left(T_{prod} - T \right) \right)\bigg\vert_{\lambda_{k-1}} \end{matrix} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] = 0 \tag{15} \]

Initial Guess for \(T_{prod}\):

\[ T_{guess} = T_{feed} + 10 \text{ K} \tag{16} \]

I also will need to write a residuals function and provide it to the ATE solver. The ATE solver will assume that it takes a single argument, a guess for the solution, and that it returns the value of the residual corresponding to that guess.

Looking at the residual expression, I see that in addition to given and known constants, it also contains the ages, \(\lambda_k\), derivatives, \(\frac{dF}{d\lambda}\bigg\vert_{\lambda_k}\), and the fluid element temperature, \(T\big\vert_{\lambda_k}\), at each of the integrand evaluation points. I will need to calculate these quantities and make them available to the residuals funtion before I call the ATE solver because they cannot be passed to the residuals function as arguments.

Residual Function:

Argument: \(T_{prod}\).

Must be Available: \(\overline{T}\), \(\overline{\lambda}\), and \(\overline{\frac{dF}{d\lambda}}\).

Returns: \(\epsilon\).

Algorithm:

\[ 0 = \epsilon = 0.5 \dot{V}_{feed} \breve{C}_p \sum_{k=2}^{N_\lambda} \left[ \left( \begin{matrix} \left( \frac{dF}{d\lambda} T_{prod} - T \right)\bigg\vert_{\lambda_k} \\+ \left( \frac{dF}{d\lambda} \left(T_{prod} - T \right) \right)\bigg\vert_{\lambda_{k-1}} \end{matrix} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] = 0 \tag{15} \]

To complete the SFR model, I need to write a computer function that solves the model equations and returns the product molar flow rates and temperature.

Before \(\dot{n}_{A,prod}\), \(\dot{n}_{Z,prod}\), and \(T_{prod}\) can be calculated, \(C_A\), \(C_Z\), \(T\), and \(\frac{dF}{d\lambda}\) must be calculated at each of the integrand evaluation points. Here, the integrand evaluation points are the ages in the tabular cumulative age distribution function.

The fluid element concentrations and temperature can be obtained by calling the fluid element function. Equations 23.8 and 23.9 can be used to calculate \(\frac{dF}{d\lambda}\) at each integrand evaluation point. At that point, the product molar flow rates of A and Z can be calculated directly.

After making \(\overline{T}\), \(\overline{\lambda}\), and \(\overline{\frac{dF}{d\lambda}}\) available to the residuals function, the product temperature can be calculated by calling an ATE solver.

SFR Function:

Arguments: none.

Must be Available: Fluid Element Function.

Returns: \(\dot{n}_{A,prod}\), \(\dot{n}_{Z,prod}\), and \(T_{prod}\).

Algorithm:

\[ \begin{matrix} \text{Fluid Element Function} \\ \Downarrow \\ \overline{C}_A, \overline{C}_Z, \overline{T} \end{matrix} \tag{17} \]

\[ \overline{\frac{dF}{d\lambda}} : \frac{dF}{d\lambda}\bigg\vert_{\lambda_k} = \frac{F_{k+1} - F_k}{\lambda_{k+1} - \lambda_k} \text{ and } \frac{dF}{d\lambda}\bigg\vert_{\lambda_f} = 0.0 \tag{18} \]

\[ \dot{n}_{A,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(C_A \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( C_A \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{13} \]

\[ \dot{n}_{Z,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(C_Z \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( C_Z \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{14} \]

\[ \overline{T}, \overline{\lambda}, \overline{\frac{dF}{d\lambda}} \, \Rightarrow \, \text{available to Resduals Function} \tag{19} \]

\[ \begin{matrix} T_{guess}, \text{Residuals Function} \\ \Downarrow \\ \text{ATE Solver} \\ \Downarrow \\ T_{prod} \end{matrix} \tag{20} \]

At this point all that remains is to calculate the quantities of interest. The narrative requests the conversion and the product temperature. The product temperature has already been calculated. To calculate the conversion,I can use its definition, Equation 3.4.

Quantities of Interest

Algorithm:

\[ f_{A,prod} = \frac{\dot{V}_{feed} C_{A,feed} - \dot{n}_{A,prod}}{\dot{V}_{feed} C_{A,feed}} \tag{21} \]

23.7.1.3 Results, Analysis and Discussion

The calculations were performed as described above. Specifically, the Derivatives Function, Fluid Element Function, Residual Function, and SFR Function were written using a mathematics software package. The product stream molar flow rates of A and Z and the product stream temperature were obtained by calling the SFR Function. The product stream conversion was then calculated using equation (21). The conversion was found to equal 93.1% and the product temperature to be 310 K. The non-ideal reactor operates adiabatically, and the reaction is exothermic. Therefore the increase in temperature is expected.

At first glance, it seems as if the 10 L volume of the non-ideal reactor was not used in the calculations, begging the question of how that is possible. The key to understanding is found in the use of the differential age distribution function to partition the feed into fluid elements as indicated in Equation 23.2. Each fluid element has a specific age, or residence time in the reactor, and it constitutes a fixed fraction of the total volumetric feed.

Dividing the volumetric feed associated with a fluid element by its residence time then yields the reactor volume associated with that fluid element. Adding up the volumes associated with all of the fluid elements (i. e. integration over all possible fluid element residence times) gives the total volume. So while the total volume does not appear explicitly in the calculations, the outlet molar flow rates and the outlet temperature are calculated by summing over the volumes of all of the fluid elements that make up the total volume.

23.7.2 An SFR Model with a Gas-Phase Reaction

The reaction between A and B, equation (1), is exothermic, \(\Delta H\) = -7.2 kcal mol-1. The rate expression is given in equation (2) with the Arrhenius pre-exponential factor, \(k_0\), equal to 1.37 x 105 m3 mol-1 min-1 and the activation energy, \(E\), to 11.1 kcal mol-1. A gas mixture containing 10% A, 65% B and 25% inert, I, flows into a non-ideal adiabatic stirred tank reactor at 165 °C and 5 atm. The feed flow rate is 1.0 m3 min-1, and the heat capacities of A, B, Z, and I are equal to 7.6, 8.2, 13.8, and 4.3 cal mol-1 K-1, respectively. The cumulative age distribution function is available in the file, cum_age_dist_fcn.csv. What conversion and temperature are predicted by an SFR model based on that age distribution?

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

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


This assignment involves a non-ideal flow reactor that will be analyzed using an SFR model. I’ll begin by summarizing the assignment. I know that I’ll need to assume a fluid element volume as a basis for the calculations, so I’ll include that in the summary. I’ll use subscripted “feed” and “prod” to indentify values for the non-ideal reactor feed and product streams, and an overbar to indicate a set of values for a quantity.

23.7.2.1 Assignment Summary

Reaction:

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

Rate Expression:

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

Reactor System: Steady-state, adiabatic, constant pressure SFR.

Quantities of Interest: \(f_A\) and \(T_{prod}\)

Given and Known Constants: \(\Delta H\) = -7.2 kcal mol-1, \(k_0\) = 1.37 x 105 m3 mol-1 min-1, \(E\) = 11.1 kcal mol-1, \(y_{A,feed}\) = 0.1, \(y_{B,feed}\) = 0.65, \(y_{I,feed}\) = 0.2, \(T_{feed}\) = 165 °C, \(P\) = 5 atm, \(\dot{V}_{feed}\) = 1.0 m3 min-1, \(\hat{C}_{p,A}\) = 7.6 cal mol-1 K-1, \(\hat{C}_{p,B}\) = 8.2 cal mol-1 K-1, \(\hat{C}_{p,Z}\) = 13.8 cal mol-1 K-1, \(\hat{C}_{p,I}\) = 4.3 cal mol-1 K-1, \(R_{en}\) = 1.987 cal mol-1 K-1, \(R_{pv}\) = 8.206 x 10-5 m3 atm mol-1 K-1, and \(\overline{F}\) vs. \(\overline{\lambda}\).

Basis: \(V_{fe}\) = 1.0 m3.

23.7.2.2 Mathematical Formulation of the Analysis

The formulation that follows assumes that the given and known constants identified in the assignment summary are available at any point in the analysis.

The “reactor” in an SFR model is a fluid element, and the BSTR design equations are used to represent it, so I’ll start by formulating the fluid element model. The reactor is adiabatic, so mole balances and an energy balance on the reacting fluid need to be included in the design equations.

The general form of the BSTR mole balance is given in Equation 6.8. Since there is only one reaction, the summation reduces to a single term, and it isn’t necessary to index the reaction.

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

The general form of the BSTR energy balance is given in Equation 6.9. Here the pressure is constant, so its time derivative is zero. The work done by the system can be assumed to be negligible, and since the reactor is adiabatic, the heat transfer term is equal to zero. The summation again reduces to a single term, and indexing the reaction isn’t necessary. Finally, it can be seen that all of the terms in the energy balance except for \(P\frac{dV}{dt}\) will have units of energy per time. The units of the \(P\frac{dV}{dt}\) term will be pressure times volume per time. The units of that term can be converted by multiplying it by the ideal gas constant in energy units, \(R_{en}\), divided by the ideal gas constant in pressure-volume units, \(R_{pv}\).

\[ \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) \]

\[ \left(\sum_i n_i \hat{C}_{p,i} \right) \frac{dT}{dt} - P\frac{R_{en}}{R_{pv}} \frac{dV}{dt} = - r V \Delta H \]

The number of dependent variables in the design equations above is one greater than the number of equations. Consequently either a dependent variable needs to be eliminated or an equation needs to be added. I’ll add the differential form of the ideal gas law, Equation 6.11 so that the number of dependent variables in the design equations equals the number of design equations. As above, the pressure is constant, so its time-derivative is zero.

\[ R_{pv} \left( T \left( \sum_i \frac{dn_i}{dt} \right) + \left( \sum_i n_i \right) \frac{dT}{dt} \right) - P \frac{dV}{dt} - V \cancelto{0}{\frac{dP}{dt}} = 0 \]

Fluid Element Model

Design Equations

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

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

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

\[ \frac{dn_I}{dt} = 0.0 \tag{6} \]

\[ \left( n_A \hat{C}_{p,A} + n_B \hat{C}_{p,B} + n_Z \hat{C}_{p,Z} + n_I \hat{C}_{p,I}\right) \frac{dT}{dt} - P\frac{R_{en}}{R_{pv}}\frac{dV}{dt} = - r V \Delta H \tag{7} \]

\[ R_{pv}T\left(\frac{dn_A}{dt} + \frac{dn_B}{dt} + \frac{dn_Z}{dt} + \frac{dn_I}{dt}\right) + R_{pv}\left( n_A + n_B + n_Z + n_I \right) \frac{dT}{dt} - P \frac{dV}{dt} = 0 \tag{8} \]

The design equations are IVODEs, and I will solve them numerically using an IVODE solver. To do that, I’ll need initial values, as stopping criterion, and a derivatives function. I can define \(t=0\) to be the instant the feed is segregated into fluid elements. That means that the initial mole fractions in the fluid element will equal the mole fractions in the non-ideal reactor feed and the initial temperature will equal the feed temperature. Having assumed the volume of the fluid element as a basis, I can use the ideal law to calculate the initial molar amounts. The feed mole fraction of Z is equal to zero, so the initial moles of Z in the fluid element is zero.

I will need to solve the design equations using each integrand evaluation point as the reaction time. Here the cumulative age distribution function is provided as a table, so the integrand evaluations points are simply the ages in the cumulative age distribution function table.

It is also important to check that the cumulative age distribution function spans the full range from \(F=0\) to \(F=1\). Examination of the table for this assignment shows that it does.

Initial Values and Stopping Criterion:

Table 23.2: Initial values and stopping criterion for solving the design equations, (3) through (8).
Variable Initial Value Stopping Criterion
\(t\) \(0\) each \(\lambda_k\) in \(\overline{\lambda}\)
\(n_A\) \(y_{A,feed}\frac{PV_{fe}}{R_{pv}T_{feed}}\)
\(n_B\) \(y_{B,feed}\frac{PV_{fe}}{R_{pv}T_{feed}}\)
\(n_Z\) \(0.0\)
\(n_I\) \(y_{I,feed}\frac{PV_{fe}}{R_{pv}T_{feed}}\)
\(T\) \(T_{feed}\)
\(V\) \(V_{fe}\)

To solve the design equations numerically, I’ll also need to write a derivatives function and provide it to the solver. Because it will be called by the solver, the only arguments that can be passed to it are the values of the independent, \(t\), and dependent variables, \(n_A\), \(n_B\), \(n_Z\), \(n_I\), \(T\), and \(V\) at the start of an integration step, and the only quantities it can return are the derivatives of the dependent variables with respect to the independent variable.

Before I can evaluate the derivatives, I need to calculate any unknown quantities appearing in the design equations. For this assignment, the only unknown in the design equations is the rate, \(r\). It can be calculated using the rate expression provided in the problem narrative. That, in turn, will require using the Arrhenius expression, Equation 4.8, to calculate the rate coefficient, and the defining equation for concentration, Equation 1.7.

I also need to write the design equations in the form of derivative expressions (see Appendix J). One way to do so is to write the design equations as a matrix equation.

\[ \boldsymbol{M} \begin{bmatrix} \frac{dn_A}{dt} \\ \frac{dn_B}{dt} \\ \frac{dn_Z}{dt} \\ \frac{dn_I}{dt} \\ \frac{dT}{dt} \\ \frac{dV}{dt} \end{bmatrix} = \underline{g} \]

Here, the mass matrix, \(\boldsymbol{M}\), and the right-hand side vector, \(\underline{g}\) are then defined as follows.

\[ \boldsymbol{M} = \begin{bmatrix} 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0\\ 0 & 0 & 0 & 0 & \displaystyle\sum_i\left(n_i \hat{C}_{p,i} \right) & -P\frac{R_{en}}{R_{pv}} \\ R_{pv}T & R_{pv}T & R_{pv}T & R_{pv}T & R_{pv}\displaystyle\sum_i n_i & -P \end{bmatrix} \]

\[ \underline{g} = \begin{bmatrix} - rV \\ - rV \\ rV \\ 0.0 \\ -rV\Delta H \\ 0 \end{bmatrix} \]

That then allows the values of the derivatives to be calculated using the inverse of the mass matrix.

\[ \begin{bmatrix} \frac{dn_A}{dt} \\ \frac{dn_B}{dt} \\ \frac{dn_Z}{dt} \\ \frac{dn_I}{dt} \\ \frac{dT}{dt} \\ \frac{dV}{dt} \end{bmatrix} = \boldsymbol{M}^{-1} \underline{g} \]

Derivatives Function:

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

Returns: \(\frac{dn_A}{dt}\), \(\frac{dn_B}{dt}\), \(\frac{dn_Z}{dt}\), \(\frac{dn_I}{dt}\), \(\frac{dT}{dt}\), and \(\frac{dV}{dt}\).

Algorithm:

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

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

\[ k = k_0 \exp{\left(\frac{-E}{RT}\right)} \tag{11} \]

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

\[ \boldsymbol{M} = \begin{bmatrix} 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0\\ 0 & 0 & 0 & 0 & \displaystyle\sum_i\left(n_i \hat{C}_{p,i} \right) & -P\frac{R_{en}}{R_{pv}} \\ R_{pv}T & R_{pv}T & R_{pv}T & R_{pv}T & R_{pv}\displaystyle\sum_i n_i & -P \end{bmatrix} \tag{12} \]

\[ \underline{g} = \begin{bmatrix} - rV \\ - rV \\ rV \\ 0.0 \\ -rV\Delta H \\ 0 \end{bmatrix}\tag{13} \]

\[ \begin{bmatrix} \frac{dn_A}{dt} \\ \frac{dn_B}{dt} \\ \frac{dn_Z}{dt} \\ \frac{dn_I}{dt} \\ \frac{dT}{dt} \\ \frac{dV}{dt} \end{bmatrix} = \boldsymbol{M}^{-1} \underline{g} \tag{14} \]

The purpose of the fluid element model is to calculate the concentration of the reagents, temperature, and relative volume of the fluid element at each of the integrand evaluation points. In this assignment, the integrand evaluation points are the ages, \(\lambda_k\), in the cumulative age distribution function.

To do that, I can set the final value of the reaction time equal to each of those ages, and solve the design equations. That will give me corresponding sets of values of \(t\), \(n_A\), \(n_B\), \(n_Z\), \(n_I\), \(T\), and \(V\) that span the range from \(t=0\) to \(t=\lambda_k\). The molar amounts and the volume at \(t=\lambda_k\) can be used to calculate the concentrations.

\[ C_i\big\vert_{t_r=\lambda_k} = \frac{n_i\big\vert_{t_r=\lambda_k}}{V\big\vert_{t_r=\lambda_k}} \]

The relative volume at \(t=\lambda_k\) is simply the volume divided by the initial volume.

\[ \gamma\big\vert_{t_r=\lambda_k} = \frac{V\vert_{t_r=\lambda_k}}{V_{fe}} \]

The reagent concentrations, fluid element temperatures, and relative volumes for each \(\lambda_k\) can then be combined into vectors with values that correspond to the ages, \(\overline{\lambda}\).

Fluid Element Function:

Arguments: none

Returns: \(\overline{C}_A\), \(\overline{C}_B\), \(\overline{C}_Z\), \(\overline{C}_I\), \(\overline{T}\), and \(\overline{\gamma}\)

Algorithm:

Repeat (15) through (12) for each \(\lambda_k\) in \(\overline{\lambda}\).

\[ \begin{matrix} \text{initial values, stopping criterion, derivatives function} \\ \Downarrow \\ \text{IVODE solver} \\ \Downarrow \\ \overline{t}, \overline{n}_A, \overline{n}_B, \overline{n}_Z, \overline{n}_I, \overline{T}, \overline{V} \end{matrix} \tag{15} \]

 

\[ C_{A,k} \text{ in } \overline{C}_A = \frac{n_A \big\vert_{t=\lambda_k}}{V\big\vert_{t=\lambda_k}} \tag{16} \]

\[ C_{B,k} \text{ in } \overline{C}_B = \frac{n_B \big\vert_{t=\lambda_k}}{V\big\vert_{t=\lambda_k}} \tag{17} \]

\[ C_{Z,k} \text{ in } \overline{C}_Z = \frac{n_Z \big\vert_{t=\lambda_k}}{V\big\vert_{t=\lambda_k}} \tag{18} \]

\[ C_{I,k} \text{ in } \overline{C}_I = \frac{n_I \big\vert_{t=\lambda_k}}{V\big\vert_{t=\lambda_k}} \tag{19} \]

\[ T_k \text{ in } \overline{T} = T\big\vert_{t=\lambda_k} \tag{20} \]

\[ \gamma_k \text{ in } \overline{\gamma} = \frac{V\big\vert_{t=\lambda_k}}{V_{fe}} \tag{21} \]

Having formulated a fluid element model, I now need to formulate a SFR model to use the results it generates. According to the SFR model, the outlet molar flow rates for a non-ideal reactor are related to the differential age distribution function and the results from the fluid element model as indicated in Equation 23.11.

\[ \dot{n}_{i,out} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(\gamma C_i \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( \gamma C_i \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \]

Before that equation can be used, I need to calculate \(\frac{dF}{d\lambda}\bigg\vert_{\lambda_k}\) at each of the integrand evaluation points. That can be done using Equations 23.8 and 23.9.

To calculate the outlet temperature for the non-ideal reactor, Equation 23.12 is used.

\[ T_{out} \,: 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left(\begin{matrix} \left(\gamma \frac{dF}{d\lambda} \left( \sum_i C_i \int_T^{T_{out}} \hat{C}_{p,i} dT \right)\right)\bigg\vert_{\lambda_k} \\+ \left(\gamma \frac{dF}{d\lambda} \left( \sum_i C_i \int_T^{T_{out}} \hat{C}_{p,i} dT \right) \right)\bigg\vert_{\lambda_{k-1}} \end{matrix} \right) \left(\lambda_k - \lambda_{k-1}\right) \right] = 0 \]

Here the heat capacities are constants so they can be factored out of the senisble heat integrals and those integrals can be evaluated analytically. Since this is an implicit equation for \(T_{out}\), it will be solved using an ATE solver, so I’ll write it as a residual expression.

\[ 0 = \epsilon = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left(\begin{matrix} \left(\gamma \frac{dF}{d\lambda} \left( \sum_i \left(C_i \hat{C}_{p,i}\right) \left(T_{out} - T \right)\right)\right)\bigg\vert_{\lambda_k} \\+ \left(\gamma \frac{dF}{d\lambda} \left( \sum_i \left(C_i \hat{C}_{p,i} \right) \left(T_{out} - T\right) \right) \right)\bigg\vert_{\lambda_{k-1}} \end{matrix} \right) \left(\lambda_k - \lambda_{k-1}\right) \right] \]

SFR Model

Model Equations

\[ \dot{n}_{A,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(\gamma C_A \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( \gamma C_A \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{22} \]

\[ \dot{n}_{B,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(\gamma C_B \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( \gamma C_B \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{23} \]

\[ \dot{n}_{Z,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(\gamma C_Z \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( \gamma C_Z \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{24} \]

\[ \dot{n}_{I,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(\gamma C_I \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( \gamma C_I \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{25} \]

\[ 0 = \epsilon = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left(\begin{matrix} \left(\gamma \frac{dF}{d\lambda} \left( \sum_i \left(C_i \hat{C}_{p,i}\right) \left(T_{prod} - T \right)\right)\right)\bigg\vert_{\lambda_k} \\+ \left(\gamma \frac{dF}{d\lambda} \left( \sum_i \left(C_i \hat{C}_{p,i} \right) \left(T_{prod} - T\right) \right) \right)\bigg\vert_{\lambda_{k-1}} \end{matrix} \right) \left(\lambda_k - \lambda_{k-1}\right) \right] \tag{26} \]

I will need to provide a guess for \(T_{prod}\) to the ATE solver. Since the reaction is exothermic and the reactor is adiabatic, I’ll guess that the outlet temperature is 10 K greater than the feed. If the solver fails to converge, I may need to adjust this guess.

I’ll also need to provide a residual function to the ATE solver. It must accept a guess for \(T_{prod}\) as its only argument and it must return the corresponding value of \(\epsilon\). Looking at the residual expression, I see that in addition to given and known constants, it also contains the derivatives, \(\frac{dF}{d\lambda}\bigg\vert_{\lambda_k}\), relative volume, \(\gamma\big\vert_{\lambda_k}\), concentrations, \(C_i\big\vert_{\lambda_k}\), and the fluid element temperature, \(T\big\vert_{\lambda_k}\), at each of the integrand evaluation points. I will need to calculate these quantities and make them available to the residuals funtion before I call the ATE solver because they cannot be passed to the residuals function as arguments.

Initial Guess for \(T_{prod}\):

\[ T_{guess} = T_{feed} + 10 \text{ K} \tag{27} \]

Residual Function:

Arguments: \(T_{prod}\)

Must be Available: \(\overline{\gamma}\), \(\overline{\frac{dF}{d\lambda}}\), \(\overline{T}\), \(\overline{C}_A\), \(\overline{C}_B\), \(\overline{C}_Z\), and \(\overline{C}_I\)

Returns: \(\epsilon\)

Algorithm:

\[ 0 = \epsilon = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left(\begin{matrix} \left(\gamma \frac{dF}{d\lambda} \left( \sum_i \left(C_i \hat{C}_{p,i}\right) \left(T_{prod} - T \right)\right)\right)\bigg\vert_{\lambda_k} \\+ \left(\gamma \frac{dF}{d\lambda} \left( \sum_i \left(C_i \hat{C}_{p,i} \right) \left(T_{prod} - T\right) \right) \right)\bigg\vert_{\lambda_{k-1}} \end{matrix} \right) \left(\lambda_k - \lambda_{k-1}\right) \right] \tag{26} \]

The last component of the SFR model is a computer function that solves the model equations and returns the product molar flow rates and temperature.

Before those quantities can be calculated, \(C_A\), \(C_B\), \(C_Z\), \(C_I\), \(T\), \(\lambda\), and \(\frac{dF}{d\lambda}\) must be calculated at each of the integrand evaluation points (i. e. the ages in the tabular cumulative age distribution function).

The fluid element concentrations, temperature, and relative volume can be obtained by calling the fluid element function. Equations 23.8 and 23.9 can be used to calculate \(\frac{dF}{d\lambda}\) at each integrand evaluation point. once that is done, the product molar flow rates can be calculated directly.

After making \(\overline{\gamma}\), \(\overline{\frac{dF}{d\lambda}}\), \(\overline{T}\), \(\overline{C}_A\), \(\overline{C}_B\), \(\overline{C}_Z\), and \(\overline{C}_I\) available to the residuals function, the product temperature can be calculated by calling an ATE solver.

SFR Function:

Arguments: none

Must be Available: Fluid Element Function

Returns: \(\dot{n}_{A,prod}\), \(\dot{n}_{B,prod}\), \(\dot{n}_{Z,prod}\), \(\dot{n}_{I,prod}\), and \(T_{prod}\).

Algorithm:

\[ \begin{matrix} \text{Fluid Element Function} \\ \Downarrow \\ \overline{C}_A, \overline{C}_B, \overline{C}_Z, \overline{C}_I, \overline{T}, and \overline{\gamma} \end{matrix} \tag{28} \]

\[ \overline{\frac{dF}{d\lambda}} : \frac{dF}{d\lambda}\bigg\vert_{\lambda_k} = \frac{F_{k+1} - F_k}{\lambda_{k+1} - \lambda_k} \text{ and } \frac{dF}{d\lambda}\bigg\vert_{\lambda_f} = 0.0 \tag{29} \]

\[ \dot{n}_{A,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(\gamma C_A \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( \gamma C_A \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{22} \]

\[ \dot{n}_{B,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(\gamma C_B \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( \gamma C_B \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{23} \]

\[ \dot{n}_{Z,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(\gamma C_Z \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( \gamma C_Z \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{24} \]

\[ \dot{n}_{I,prod} = 0.5 \dot{V}_{feed} \sum_{k=2}^{N_\lambda} \left[ \left( \left(\gamma C_I \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_k} + \left( \gamma C_I \frac{dF}{d\lambda}\right)\bigg\vert_{\lambda_{k-1}} \right)\left(\lambda_k - \lambda_{k-1}\right) \right] \tag{25} \]

\[ \overline{\gamma}, \overline{\frac{dF}{d\lambda}}, \overline{T}, \overline{C}_A, \overline{C}_B, \overline{C}_Z, \overline{C}_I \, \Rightarrow \, \text{available to Residuals Function} \tag{30} \]

\[ \begin{matrix} T_{guess}, \text{Residuals Function} \\ \Downarrow \\ \text{ATE Solver} \\ \Downarrow \\ T_{prod} \end{matrix} \tag{31} \]

At this point all that remains is to calculate the quantities of interest. The narrative requests the conversion and the product temperature. The product temperature has already been calculated. I can use the definition of conversion, Equation 3.4, to calculate its value.

Quantities of Interest

Algorithm:

\[ \dot{n}_{A,feed} = y_{A,feed} \frac{P \dot{V}_{feed}}{R_{pv}T_{feed}} \tag{32} \]

\[ f_{A,prod} = \frac{\dot{n}_{A,feed} - \dot{n}_{A,prod}}{\dot{n}_{A,feed}} \tag{33} \]

23.7.2.3 Results, Analysis and Discussion

The calculations were performed as described above. Specifically, the Derivatives Function, Fluid Element Function, Residual Function, and SFR Function were written using a mathematics software package. The product stream molar flow rates of A, B, Z, and I and the product stream temperature were obtained by calling the SFR Function. The product stream conversion was then calculated using equations (32) and (33). The conversion was found to equal 78.4% and the product temperature to be 191 °C.

The temperature increase from 165 to 191 °C is expected. The reaction is exothermic, so heat is released. The reactor is adiabatic, so the released heat is not removed and instead causes the gas temperature to rise.

23.8 Symbols Used in Chapter 23

Symbol Meaning
\(dF\) Differential age distribution function.
\(f_i\) Conversion of reagent \(i\).
\(\overline{g}\) Right-hand side vector when ODEs are written as a matrix equation.
\(k\) Rate coefficient.
\(k_0\) Arrhenius pre-exponential factor.
\(n_i\) Molar amount of reagent \(i\); a subscripted “s,0” denotes the initial sample and “s,f” denotes the sample after batch processing.
\(\dot{n}_i\) Molar flow rate; subscripted “feed” and “product” denote the flow stream, a subscripted “\(\tau\)” denotes a reaction time of \(\tau\).
\(r_j\) Rate of reaction \(j\).
\(t\) Time, a subscripted “r” denotes reaction time.
\(y\) Integrand or mole fraction.
\(C_i\) Concentration of reagent \(i\).
\(\hat{C}_{p,i}\) Molar heat capacity of reagent \(i\).
\(\breve{C}_p\) Volumetric heat capacity.
\(E\) Activation energy.
\(F\left(\lambda\right)\) Cumulative age distribution function.
\(\boldsymbol{M}\) Mass matrix.
\(P\) Pressure.
\(R\) Ideal gas constant.
\(T\) Temperature.
\(V\) Reacting fluid volume.
\(\dot{V}\) Volumetric flow rate; subscripted “feed” and “product” denote the flow stream, a subscripted “\(\tau\)” denotes a reaction time of \(\tau\).
\(\epsilon\) Residual.
\(\gamma\) Relative fluid volume (current volume divided by initial volume).
\(\lambda\) Fluid age.
\(\rho\) Density.
\(\tau\) Reaction time.
\(\Delta H\) Heat of reaction.
\(\Delta t_B\) Sampling time used as a basis for calculations.