# Optimization of the Solution of a Dispersion Model

^{1}

^{2}

^{3}

^{4}

^{5}

^{*}

Previous Article in Journal

Department of Hydraulics and Environment Protection, Technical University of Civil Engineering of Bucharest, Bd. Lacul Tei 122-124, 020396 Bucharest, Romania

Department of Mathematics and Computer Science, Technical University of Civil Engineering of Bucharest, Bd. Lacul Tei 122-124, 020396 Bucharest, Romania

Department of Mechanical Technology, Technical University of Civil Engineering of Bucharest, Bd. Lacul Tei 122-124, 020396 Bucharest, Romania

Department of Electrical Engineering, Automation and Informatics, Faculty of Engineering, Slovak University of Agriculture in Nitra, Tr. A. Hlinku 2, 949 76 Nitra, Slovakia

Department of Mechanical Engineering, Faculty of Technology, Institute of Technology and Business in České Budějovice, Okružní 10, 370 01 České Budějovice, Czech Republic

Author to whom correspondence should be addressed.

Received: 3 February 2020 / Revised: 19 February 2020 / Accepted: 21 February 2020 / Published: 1 March 2020

The study of the combination of chemical kinetics with transport phenomena is the main step for reactor design. It is possible to deviate from the model behaviour, the cause of which may be fluid channelling, fluid recirculation, or creation of stagnant regions in the vessel, by using a dispersion model. In this paper, the known general solution of the dispersion model for closed vessels is given in a new, straightforward form. In order to improve the classical theoretical solution, a hybrid of analytical and numerical methods is used. It is based on the general analytic solution and the least-squares method by fitting the results of a tracer test carried out on the vessel with the values of the analytic solution. Thus, the accuracy of the estimation for the vessel dispersion number is increased. The presented method can be used to similar problems modelled by a partial differential equation and some boundary conditions which are not sufficient to ensure the uniqueness of the solution.

The efficiency of ensuring stable and economically optimal operating regimes for a wastewater treatment plant requires a mathematical model that describes, using specific tools, the simulation of the technological processes that occur in its operation. As the main technological object of a treatment plant is represented by the biological reactor, the modelling of the complex phenomena that happen in the wastewater flow process, with the help of mathematical models, is a real aid in understanding its functioning, as well as in optimizing the design and the energy consumption.

The study of the combination of chemical kinetics with transport phenomena is the main step for reactor design. Inside of a reactor, the flow pattern usually deviates from the ideal cases of plug flow or mixed flow. This deviation, the cause of which may be fluid channelling, recirculation of fluid, or the fact that stagnant regions were created in the vessel, can be approached by using a model characterized by one parameter (e.g., tanks-in-series model or dispersion model) or a model characterized by two parameters (e.g., reactor with bypassing and dead volume). The parameter(s) of the model is (are) evaluated using the residence time distribution (RTD). MacMullin and Weber were the first who came up with the idea of using the RTD [1]. In the 1950s, Danckwerts proposed the concept of RTD [2], which has been used by many authors. It is mostly associated with the analysis of the chemical reactor performance [3,4,5,6]. A careful analysis of the RTD theory has been published by Nauman [7], who introduced new modelling techniques and applications. The concept of RTD has some limitations mainly for complex systems. In the process simulation of single flow pattern in wastewater treatment, the main restriction is that the process is running under steady-state conditions during the tracer experiment [8]. Therefore, some attempts to use this approach in the transient state were made [9,10,11]. For more than 70 years, this method has been used for process optimization [12,13,14,15,16,17,18].

When we introduce a model pulse of tracer into the fluid that enters a vessel, the pulse spreads while passing through the vessel. If we want to characterize the spreading based on this model, a diffusion-like process superimposed on plug flow, called dispersion or longitudinal dispersion, is assumed [4]. This spreading process is represented by the dispersion coefficient D (${m}^{2}\cdot {s}^{-1}$). A high value of D denotes the tracer curve rapid spreading, a low value of D represents slow spreading, and D = 0 implies that there is no spreading, hence plug flow. Moreover, the vessel dispersion number, which is a form of the inverse of the Peclet number Pe^{−1}, D/µL where u is the mean (displacement) velocity and L is the length of the vessel, is the dimensionless group characterizing the spread in the whole vessel. If back mixing is transferred on a plug flow of a fluid to a certain degree, whose magnitude does not depend on the position within the vessel, it is called the dispersion model [4]. To solve the problem of dispersion, minimization is a challenging task. There is an extensive literature on the dispersion model [19,20,21,22]. In fact, the dispersion model is based on the dispersion Equation (1) and two boundary conditions, Equations (2) and (3). This problem has infinitely many solutions, but by using additional ideal theoretical conditions, Murphy and Timpany [23] selected an analytic solution described by the residence time distribution (see Equation (33)). However, in the real cases, some errors may often appear by comparing the analytic solution, Equation (33), with the results obtained by measurements (e.g., by so-called tracer test).

In the paper, we introduce a new method for estimating the vessel dispersion number. The outline of this paper is as follows: In Section 2.1, the general solution of the dispersion model for closed vessels is obtained by separation of variables to underline the possibility to improve the known solution obtained by Murphy and Timpany in Reference [23]. Here, a simplified expression of the general solution and of the residence time distribution E(t) is also obtained. Section 2.2 deals with a new solution which is achieved by fitting the results of a tracer test carried out on the vessel with the values of the analytic solution. In Section 2.3, an iterative method to improve the accuracy of the estimation of the parameter D/µL is presented. In fact, an inverse problem is solved by using the results of a tracer test and the least-squares method. An illustrative numerical example is given, and the conclusions of the work are presented in Section 3. The standard tools used for partial differential equations can be found, for example, in Reference [24].

Finally, the given approach is an excellent example of interdisciplinary mathematics with a direct impact on the needs of current industrial practice.

The basic differential equation representing the dispersion model is (see Reference [23,25] or Reference [4]),
where $c(x,t)$ is the concentration (wt./unit vol.) at time t and distance x. We have to take two cases into consideration: boundary conditions for closed vessels and open vessels called Danckwerts boundary conditions (see Reference [3], p. 883). In closed vessels, no dispersion or radial variation in concentration is assumed upstream (closed) nor downstream (closed) of the section where the reaction occurs; therefore, this is a closed-closed vessel. In open vessels, dispersion occurs both upstream (open) and downstream (open) of the reaction section; therefore, this is an open-open vessel. In this paper, we discuss the first case that models all the reactors for which there is no flow or diffusion or up-flow eddies at the entrance or at the vessel exit (see Reference [4]).

$$D\frac{{\partial}^{2}c}{\partial {x}^{2}}-u\frac{\partial c}{\partial x}=\frac{\partial c}{\partial t},$$

Thus, in the case of a closed vessel, plug flow (no dispersion) can be found to the immediate left of the entrance and the immediate right of the exit, so that the boundary conditions can be written in the following form (see Reference [25] or Reference [3], p. 884):

$$\frac{\partial c}{\partial x}(0,t)=\frac{u}{D}c(0,t),$$

$$\frac{\partial c}{\partial x}(L,t)=0.$$

In order to find a solution of Equation (1) that verifies the boundary conditions, Equations (2) and (3) (i.e., to solve the boundary value problem (1)–(3)), the method of separation of variables is used. Thus, by choosing
from Equation (1) it follows that
where $\lambda $ is a constant. Hence, the following differential equations are obtained:
and

$$c(x,t)=X(x)T(t)$$

$$\frac{{X}^{\u2033}(x)-\frac{u}{D}{X}^{\prime}(x)}{X(x)}=\frac{{T}^{\prime}(t)}{DT(t)}=\lambda ,$$

$${T}^{\prime}(t)-\lambda DT(t)=0$$

$${X}^{\u2033}(x)-\frac{u}{D}{X}^{\prime}(x)-\lambda X(x)=0$$

The general solution of the differential Equation (6) is
where ${\alpha}_{1}$ is a constant. Since $\underset{t\to \infty}{\mathrm{lim}}c(x,t)=0,$ by (8), it follows that $\lambda <0$.

$$T(t)={\alpha}_{1}{e}^{\lambda Dt},$$

By the Equation (4) and using the boundary conditions, Equations (2) and (3), it follows that
and

$${X}^{\prime}(L)=0$$

$${X}^{\prime}(0)=\frac{u}{D}X(0)$$

It is easy to prove that the conditions of Equations (9) and (10) can be fulfilled only in the case when ${\left(\frac{u}{D}\right)}^{2}+4\lambda <0$. If we denote by ${\left(\frac{u}{D}\right)}^{2}+4\lambda =-4{\delta}^{2}$, where $\delta $ is a positive number, then the general solution of the differential Equation (7) is:
where ${\alpha}_{2}$ and ${\alpha}_{3}$ are constants. Hence, the derivative ${X}^{\prime}(x)$, is written:

$$X(x)={e}^{\frac{ux}{2D}}\left({\alpha}_{2}\mathrm{cos}\delta x+{\alpha}_{3}\mathrm{sin}\delta x\right)$$

$${X}^{\prime}(x)={e}^{\frac{ux}{2D}}\cdot \frac{u}{2D}\left({\alpha}_{2}\mathrm{cos}\delta x+{\alpha}_{3}\mathrm{sin}\delta x\right)+{e}^{\frac{ux}{2D}}\cdot \delta \left(-{\alpha}_{2}\mathrm{sin}\delta x+{\alpha}_{3}\mathrm{cos}\delta x\right)$$

Consequently, from (9)–(12), the following system is obtained:

$$\{\begin{array}{c}\left(\frac{u}{D}\mathrm{cos}\delta L-2\delta \mathrm{sin}\delta L\right){\alpha}_{2}+\left(\frac{u}{D}\mathrm{sin}\delta L+2\delta \mathrm{cos}\delta L\right){\alpha}_{3}=0\\ \frac{u}{D}{\alpha}_{2}-2\delta {\alpha}_{3}=0\end{array}$$

From the last equation of the system, Equation (13), it follows that
and, by replacing in the first one, we get:

$${\alpha}_{2}=\frac{2D\delta}{u}{\alpha}_{3}$$

$$\left(2\mathsf{\delta}\mathrm{cos}\mathsf{\delta}L-\frac{4{\mathsf{\delta}}^{2}D}{u}\mathrm{sin}\mathsf{\delta}L+\frac{u}{D}\mathrm{sin}\mathsf{\delta}L+2\mathsf{\delta}\mathrm{cos}\mathsf{\delta}L\right){\alpha}_{3}\text{}=\text{}0$$

Thus, in order to obtain nonzero solutions of the system, Equation (13), the following condition must be verified:

$$\left(\frac{2u}{D}+4\frac{\lambda D}{u}\right)\mathrm{sin}\delta L+4\delta \mathrm{cos}\delta L=0.$$

This condition may be written in the following form:

$$\mathrm{cot}\mathsf{\delta}L=-\frac{\frac{u}{D}+\frac{2\lambda D}{u}}{\sqrt{-\frac{{u}^{2}}{{D}^{2}}-4\lambda}}.$$

By denoting

$$\mu =\frac{L}{2}\sqrt{-\frac{{u}^{2}}{{D}^{2}}-4\lambda}$$

Equation (17) becomes

$$\mathrm{cot}\mu =\frac{1}{2}\left(\frac{2D\mu}{uL}-\frac{uL}{2D\mu}\right).$$

This can be written in a simplified form if we denote by $U=\frac{Lu}{2D}$,

$$\mathrm{cot}\mu =\frac{1}{2}\left(\frac{\mu}{U}-\frac{U}{\mu}\right).$$

As can be noticed in Figure 1, Equation (20) has an infinity of solutions ${\mu}_{1},\text{}{\mu}_{2},\text{}\dots ,\text{}{\mu}_{n},\text{}\dots $ which are the abscissas of the intersection points between the graphs of the functions ${f}_{1}(\mu )=\mathrm{cot}\mu $ and ${f}_{2}(\mu )=\frac{1}{2}\left(\frac{\mu}{U}-\frac{U}{\mu}\right),$ where $\mu \in \mathbf{R}\backslash \{k\pi :k\in \mathbf{Z}\}$. As a matter of fact, there is exactly one solution ${\mu}_{k}$ on each interval $((k-1)\pi ,k\pi )$, $k\in \mathbf{Z}.$

Now, by Equations (14) and (18), it follows that
and, by Equation (9), we obtain
or, equivalently,

$${\alpha}_{2}=\frac{2{\mu}_{n}D}{uL}{\alpha}_{3}=\frac{{\mu}_{n}}{U}{\alpha}_{3},$$

$${X}_{n}(x)={e}^{\frac{ux}{2D}}\left(\frac{2{\mu}_{n}D}{uL}\mathrm{cos}\frac{{\mu}_{n}x}{L}+\mathrm{sin}\frac{{\mu}_{n}x}{L}\right){\alpha}_{3}$$

$${X}_{n}(x)=\frac{2D}{uL}{e}^{\frac{ux}{2D}}\left({\mu}_{n}\mathrm{cos}\frac{{\mu}_{n}x}{L}+\frac{uL}{2D}\mathrm{sin}\frac{{\mu}_{n}x}{L}\right){\alpha}_{3}.$$

From Equations (4), (8) and (23), the following solutions of the boundary problem, Equations (1)–(3), are obtained:
where ${\gamma}_{n}$, with $n=1,2,\dots $, are constants. By denoting

$${c}_{n}(x,t)={\gamma}_{n}{e}^{\frac{ux}{2D}-\left(\frac{{u}^{2}}{4{D}^{2}}+\frac{{\mu}_{n}^{2}}{{L}^{2}}\right)Dt}\left({\mu}_{n}\mathrm{cos}\frac{{\mu}_{n}x}{L}+\frac{uL}{2D}\mathrm{sin}\frac{{\mu}_{n}x}{L}\right),$$

$${k}_{n}^{2}=\frac{{u}^{2}}{4D}+\frac{{\mu}_{n}^{2}D}{{L}^{2}},$$

Equation (24) becomes

$${c}_{n}(x,t)={\gamma}_{n}{e}^{\frac{ux}{2D}-{k}_{n}^{2}t}\left({\mu}_{n}\mathrm{cos}\frac{{\mu}_{n}x}{L}+\frac{uL}{2D}\mathrm{sin}\frac{{\mu}_{n}x}{L}\right).$$

Hence, every finite sum or a suitable infinite series whose terms have the form, Equation (26), is a solution of the boundary value problem, Equations (1)–(3). Thus, the general solution is
where m is either a positive integer or $m=+\infty .$ In the second case, also discussed in Reference [23], the constants ${\gamma}_{n}$ are small enough such that the series and its partial derivatives converge uniformly on $[0,L]\times [0,+\infty ).$

$$c(x,t)={\displaystyle \sum _{n=1}^{m}{\gamma}_{n}{e}^{\frac{ux}{2D}-{k}_{n}^{2}t}\left({\mu}_{n}\mathrm{cos}\frac{{\mu}_{n}x}{L}+\frac{uL}{2D}\mathrm{sin}\frac{{\mu}_{n}x}{L}\right)},$$

Let $\overline{t}$ be the mean residence time, which can be written (for a closed system) as $\overline{t}=\frac{L}{u}$. Since

$$\frac{L}{u}{k}_{n}^{2}=\frac{L}{u}\left(\frac{{u}^{2}}{4D}+\frac{{\mu}_{n}^{2}D}{{L}^{2}}\right)=\frac{U}{2}+\frac{{\mu}_{n}^{2}}{2U}$$

Equation (27) can be written in the equivalent form:

$$c(L,t/\overline{t})={\displaystyle \sum _{n=1}^{m}{\gamma}_{n}{e}^{U-\frac{{U}^{2}+{\mu}_{n}^{2}}{2U}\cdot \frac{t}{\overline{t}}}\left(U\mathrm{sin}{\mu}_{n}+{\mu}_{n}\mathrm{cos}{\mu}_{n}\right)}.$$

By choosing ${\gamma}_{n}={\beta}_{n}C$, where

$${\beta}_{n}:=\frac{2{\mu}_{n}}{{U}^{2}+2U+{\mu}_{n}^{2}},$$

$$T=\frac{t}{\overline{t}},\text{}\mathrm{and}$$

$$C={\displaystyle \underset{0}{\overset{\infty}{\int}}c(L,T)dT}$$

It follows Equation (10) from Reference [23] expressing the (normalized) residence time distribution $E(T)$:

$$E(T)=c(L,T)/C=2{\displaystyle \sum _{n=1}^{\infty}\frac{{\mu}_{n}(U\mathrm{sin}{\mu}_{n}+{\mu}_{n}\mathrm{cos}{\mu}_{n})}{{U}^{2}+2U+{\mu}_{n}^{2}}{e}^{U-\frac{{U}^{2}+{\mu}_{n}^{2}}{2U}\cdot T}}$$

By taking into account that, for $n=1,2,\dots $, ${\mu}_{n}$ is the unique solution of Equation (20) on the interval $((n-1)\pi ,n\pi ),$, we obtain that

$$\mathrm{sin}{\mu}_{n}:={(-1)}^{n-1}\frac{2U{\mu}_{n}}{{\mu}_{n}^{2}+{U}^{2}},\phantom{\rule{1.7cm}{0ex}}\mathrm{cos}{\mu}_{n}:={(-1)}^{n-1}\frac{{\mu}_{n}^{2}-{U}^{2}}{{\mu}_{n}^{2}+{U}^{2}}.$$

It follows that

$$U\mathrm{sin}{\mu}_{n}+{\mu}_{n}\mathrm{cos}{\mu}_{n}={(-1)}^{n-1}{\mu}_{n}$$

Hence, we obtain a new, simplified expression for the general solution $c(L,T)$ and for the (normalized) residence time distribution $E(T)$. Equations (29) and (33) can be written in the following form:
and, respectively,

$$c(L,T)={\displaystyle \sum _{n=1}^{m}{(-1)}^{n-1}{\gamma}_{n}{\mu}_{n}{e}^{U-\frac{{U}^{2}+{\mu}_{n}^{2}}{2U}\cdot T}}$$

$$E(T)=2{\displaystyle \sum _{n=1}^{\infty}\frac{{(-1)}^{n-1}{\mu}_{n}^{2}}{{U}^{2}+2U+{\mu}_{n}^{2}}{e}^{U-\frac{{U}^{2}+{\mu}_{n}^{2}}{2U}\cdot T}}$$

In this subsection, we aim to find a solution of the boundary value problem, Equations (1)–(3), more accurate than Equation (36). Consider a sequence of r pairs of real values for $(T,c(L,{T}_{i}))$:
where $c(L,T)$ is the function defined by Equation (36), ${\gamma}_{n}$ are adjustable parameters, and m is a suitable positive integer. The values of ${\gamma}_{1},{\gamma}_{2},\dots ,{\gamma}_{m}$ will be obtained by fitting the pairs (38) with $({T}_{i},c(L,{T}_{i}))$. Thus, a new solution $\overline{c}(x,t)$ will be found by considering

$$({T}_{1},\text{}{c}_{1}),\text{}({T}_{2},\text{}{c}_{2}),\text{}\dots ,({T}_{r},\text{}{c}_{r})$$

$$\overline{c}(L,T)={\displaystyle \sum _{n=1}^{m}{(-1)}^{n-1}{\gamma}_{n}{\mu}_{n}{e}^{U-\frac{{U}^{2}+{\mu}_{n}^{2}}{2U}\cdot T}}$$

As usual in the least-squares method [26], the quantity
must be minimized. Therefore, by using Equation (39), the partial derivatives of the function
with respect to ${\gamma}_{1},\dots ,{\gamma}_{m}$, vanish. Since
the following system of equations is obtained:
where $j=1,2,\dots ,m$. By solving the system, Equation (43), the constants ${\gamma}_{1},\dots ,{\gamma}_{m}$ are obtained and then, by replacing in Equation (39), the solution $\overline{c}(L,T)$ is determined.

$${\chi}^{2}={\displaystyle \sum _{i=1}^{r}{\left(\overline{c}(L,{T}_{i})-{c}_{i}\right)}^{2}}$$

$${\chi}^{2}\left({\gamma}_{1},{\gamma}_{2},\dots ,{\gamma}_{m}\right)={\displaystyle \sum _{i=1}^{r}{\left({\displaystyle \sum _{n=1}^{m}{(-1)}^{n-1}{\gamma}_{n}{\mu}_{n}{e}^{U-\frac{{U}^{2}+{\mu}_{n}^{2}}{2U}\cdot {T}_{i}}}-{c}_{i}\right)}^{2}}$$

$$\frac{\partial {\chi}^{2}}{\partial {\gamma}_{j}}=2{\displaystyle \sum _{i=1}^{r}\left(\left({\displaystyle \sum _{n=1}^{m}{\gamma}_{n}{(-1)}^{n-1}{\mu}_{n}{e}^{U-\frac{{U}^{2}+{\mu}_{n}^{2}}{2U}\cdot {T}_{i}}}-{c}_{i}\right){(-1)}^{j-1}{\mu}_{j}{e}^{U-\frac{{U}^{2}+{\mu}_{j}^{2}}{2U}\cdot {T}_{i}}\right)},\text{}j=1,2,\dots ,m,$$

$$\sum _{n=1}^{m}{(-1)}^{n-1}{\mu}_{n}\left({\displaystyle \sum _{i=1}^{r}{e}^{U-\frac{2{U}^{2}+{\mu}_{n}^{2}+{\mu}_{j}^{2}}{2U}\cdot {T}_{i}}}\right)}\text{}{\gamma}_{n}={\displaystyle \sum _{i=1}^{r}{c}_{i}{e}^{-\frac{{U}^{2}+{\mu}_{j}^{2}}{2U}\cdot {T}_{i}}},$$

In the case of a closed vessel, the following equation holds as in Reference [3,9,27], or Reference [4], p. 300:
where ${\sigma}_{t}^{2}$ is the variance. Hence, it follows an estimation of the parameter $D/uL$. In this section we present a new iterative method to increase the accuracy of this estimation. An initial value ${\left(\frac{D}{uL}\right)}_{0}$, determined by Equation (44) using the sequence in Equation (38), is considered. Let m be a fixed positive number less than r. Then, by the system in Equation (43), the values ${\gamma}_{1}^{(0)},\dots ,{\gamma}_{m}^{(0)}$ are obtained, and by replacing in Equation (39), the solution ${\overline{c}}_{0}(L,T)$ is obtained. The main idea is to find a value of ${\left(\frac{D}{uL}\right)}_{\ast}$ such that the corresponding solution ${\overline{c}}_{\ast}(L,T)$ provides the minimum value for the function

$$\frac{{\sigma}_{t}^{2}}{{\overline{t}}^{2}}=2\frac{D}{uL}-2{\left(\frac{D}{uL}\right)}^{2}\left(1-\mathrm{exp}\left(-\frac{uL}{D}\right)\right),$$

$$E{r}_{\overline{c}}\left(\frac{D}{uL}\right)={\displaystyle \sum _{i=1}^{r}{\left(\overline{c}(L,{T}_{i})-{c}_{i}\right)}^{2}}.$$

Since the expression of this function obtained by Equations (39) and (43) is complicated enough to find ${\left(\frac{D}{uL}\right)}_{\ast}$ by analytic methods, a suitable positive small real number $h$ is considered. Then, ${\left(\frac{D}{uL}\right)}_{1}$ is chosen to be that number from the set
for which the function $E{r}_{\overline{C}}\left(\frac{D}{uL}\right)$ has the minimum value. If ${\left(\frac{D}{uL}\right)}_{1}={\left(\frac{D}{uL}\right)}_{0}$, then ${\left(\frac{D}{uL}\right)}_{\ast}$ can be approximated by ${\left(\frac{D}{uL}\right)}_{0}$. In the case when ${\left(\frac{D}{uL}\right)}_{1}={\left(\frac{D}{uL}\right)}_{0}-h$, the following sequences of values is considered:

$$\left\{{\left(\frac{D}{uL}\right)}_{0}-h,\text{}{\left(\frac{D}{uL}\right)}_{0},\text{}{\left(\frac{D}{uL}\right)}_{0}+h\right\}$$

$${\left(\frac{D}{uL}\right)}_{s}={\left(\frac{D}{uL}\right)}_{0}-sh,s=2,3,\dots $$

Then, ${\left(\frac{D}{uL}\right)}_{\ast}$ is approximated by the first term ${\left(\frac{D}{uL}\right)}_{s}$, for which

$$E{r}_{\overline{C}}{\left(\frac{D}{uL}\right)}_{s+1}>E{r}_{\overline{C}}{\left(\frac{D}{uL}\right)}_{s}$$

Similarly, if ${\left(\frac{D}{uL}\right)}_{1}={\left(\frac{D}{uL}\right)}_{0}+h$, then ${\left(\frac{D}{uL}\right)}_{\ast}$ is approximated by the first term ${\left(\frac{D}{uL}\right)}_{s}={\left(\frac{D}{uL}\right)}_{0}+sh$ such that

$$E{r}_{\overline{C}}\left({\left(\frac{D}{uL}\right)}_{0}+(s+1)h\right)>E{r}_{\overline{C}}\left({\left(\frac{D}{uL}\right)}_{0}+sh\right).$$

A first-order reaction takes place in a tubular reactor 6.36 m in length with a diameter of 10 cm (see Reference [3], (p. 889)). The value of the specific reaction rate is 0.25 min^{−1}. Table 1 shows the results of a tracer test performed on this reactor.

By considering a dispersion model of a closed vessel, a comparison between the known analytic method and the new iterative method is given.

In this case $r=13$ and by Equation (44) it follows that $D/uL=1/7.5\approx 0.13$, and $U=3.75$. Since $\overline{t}=\frac{{\displaystyle \sum _{i=2}^{13}{t}_{i}{c}_{i}({t}_{i+1}-{t}_{i})}}{{\displaystyle \sum _{i=2}^{13}{c}_{i}({t}_{i+1}-{t}_{i})}}=5.0541$ (see Reference [4], p. 294), we obtain for ${T}_{i}$ the values listed in Table 2.

Then, for $m=4$, by Equations (20), (30) and (43), we obtain for ${\mu}_{i}$, ${\beta}_{i}$ and ${\gamma}_{i}$, $i=1,\dots ,4$ the values given in Table 3.

The coefficients ${\gamma}_{i}$, i = 1, …, 4 from Table 3 are used in Equation (39) to obtain the solution $\overline{c}(L,T)$. The constant C is calculated by using $\overline{c}(L,T)$ in Equation (32). The value obtained is $C=9.2099$. New coefficients ${\gamma}_{n}$ are calculated afterwards, by the formula ${\gamma}_{n}=\mathrm{C}{\beta}_{n}$. By replacing them in Equation (33), the function $c(L,T)=C\cdot E(T)$ is obtained. By comparing (using Equation (45)) the errors obtained for each solution we get that $E{r}_{c}\left(0.13\right)={\displaystyle \sum _{i=1}^{13}{\left(c(L,{T}_{i})-{c}_{i}\right)}^{2}}=33.9486$ and $E{r}_{\overline{c}}\left(0.13\right)={\displaystyle \sum _{i=1}^{13}{\left(\overline{c}(L,{T}_{i})-{c}_{i}\right)}^{2}}=4.4090$. Hence, the new solution $\overline{c}(L,T)$ leads to a better accuracy than $c(L,T).$ In order to improve the accuracy of the new theoretical model, the initial value ${\left(\frac{D}{uL}\right)}_{0}=0.13$ is considered. By choosing $h=0.01$, and ${\left(\frac{D}{uL}\right)}_{-h}=0.12$, ${\left(\frac{D}{uL}\right)}_{+h}=0.14$, and using the method from Section 2.2, the corresponding solutions ${\overline{c}}_{-h}(L,T)$ and ${\overline{c}}_{+h}(L,T)$ are determined. Then, $E{r}_{{\overline{c}}_{-h}}\left(0.12\right)=5.5915$ and $E{r}_{{\overline{c}}_{+h}}\left(0.14\right)=3.9981$ are obtained. Hence, we consider ${\left(\frac{D}{uL}\right)}_{1}=0.14,$ ${\left(\frac{D}{uL}\right)}_{2}=0.15,$ ${\left(\frac{D}{uL}\right)}_{3}=0.16,$…. and calculate the values of the error function for each case. These values are given in Table 4.

The minimum error is obtained for ${\left(\frac{D}{uL}\right)}_{\ast}\approx {\left(\frac{D}{uL}\right)}_{7}=0.20.$ In this case $E{r}_{{c}_{7}}\left(0.20\right)=9.5359$ and the graphs of ${\overline{c}}_{7}(L,T)$ and ${c}_{7}(L,T)$ are represented in Figure 2.

A comparison of the errors between the real values and those obtained by the classical solution and the new solution is presented in Figure 3.

Thus, by this new method, the solution given by Murphy and Timpany [23] can be improved by matching the results of a tracer test performed on the vessel with the values of the analytic solution. Moreover, as a solution of an inverse problem, by an iterative method, the accuracy of the estimation of the vessel dispersion number increased. The presented method can be used to similar problems modelled by a partial differential equation when the known boundary conditions are not sufficient to ensure the uniqueness of the solution.

Conceptualization, A.-N.D. and G.G.; methodology, M.J.; validation, M.H., and J.V.; formal analysis, L.R.; investigation, A.-N.D. and L.R.; data curation, S.P.; writing—original draft preparation, A.-N.D.; writing—review and editing, G.G., S.P., M.H. and J.V.; visualization, M.J.; supervision, L.R. All authors have read and agreed to the published version of the manuscript.

This research received no external funding.

The second author would like to dedicate this paper to the memory of Eng. Alexandru Braha.

The authors declare no conflict of interest.

- Mac Mullin, R.B. The Theory of Short Circuiting in Continuous-Flow Mixing Vessels in Series and the Kinetics of Chemical Reactions in Such Systems. Trans. Am. Inst. Chem. Eng.
**1935**, 31, 409–458. [Google Scholar] - Danckwerts, P.V. Continuous flow systems: Distribution of residence times. Chem. Eng. Sci.
**1953**, 2, 1–13. [Google Scholar] [CrossRef] - Fogler, H.S. Elements of Chemical Reaction Engineering, 5th ed.; Prentice Hall: Boston, MA, USA, 2016. [Google Scholar]
- Levenspiel, O. Chemical Reaction Engineering, 3rd ed.; Wiley: New York, NY, USA, 1999. [Google Scholar]
- Nauman, E.B. Chemical Reactor Design, Optimization, and Scaleup, 2nd ed.; Wiley: Hoboken, NJ, USA, 2008. [Google Scholar]
- Froment, G.F.; Bischoff, K.B.; De Wilde, J. Chemical Reactor Analysis and Design, 3rd ed.; Wiley: Hoboken, NJ, USA, 2011. [Google Scholar]
- Nauman, E.B. Residence Time Theory. Ind. Eng. Chem. Res.
**2008**, 47, 3752–3766. [Google Scholar] [CrossRef] - Matko, T.; Fawcett, N.; Sharp, A.; Stephenson, T. A numerical model of flow in circular sedimentation tanks. Process Saf. Environ.
**1996**, 74, 197–204. [Google Scholar] [CrossRef] - Nauman, E.B. Residence time distribution theory for unsteady stirred tank reactors. Chem. Eng. Sci.
**1969**, 24, 1461–1470. [Google Scholar] [CrossRef] - Niemi, A.J.; Zenger, K.; Thereska, J.; Martinez, J.G. Tracer testing of processes under variable flow and volume. Nukleonika
**1998**, 43, 73–94. [Google Scholar] - Fernandez-Sempere, J.; Font-Montesinos, R.; Espejo-Alcaraz, O. Residence time distribution for unsteady-state systems. Chem. Eng. Sci.
**1995**, 50, 223–230. [Google Scholar] [CrossRef] - Claudel, S.; Leclerc, J.P.; Tétar, L.; Lintz, H.G.; Bernard, A. Recent extensions of the residence time distribution concept: Unsteady state conditions and hydrodynamic model developments. Braz. J. Chem. Eng.
**2000**, 17, 947–954. [Google Scholar] [CrossRef] - Bogatykh, I.; Osterland, T. Characterization of Residence Time Distribution in a Plug Flow Reactor. ChemIngTech
**2019**, 91, 668–672. [Google Scholar] [CrossRef] - Hohmann, L.; Schmalenberg, M.; Prasanna, M.; Matuschek, M.; Kockmann, N. Suspension flow behavior and particle residence time distribution in helical tube devices. Chem. Eng. J.
**2019**, 360, 1371–1389. [Google Scholar] [CrossRef] - Chen, K.; Bachmann, P.; Bück, A.; Jacob, M.; Tsotsas, E. CFD simulation of particle residence time distribution in industrial scale horizontal fluidized bed. J. Powder Technol.
**2019**, 345, 129–139. [Google Scholar] [CrossRef] - Deshmukh, S.S.; Sathe, M.J.; Joshi, J.B.; Koganti, S.B. Residence time distribution and flow patterns in the single-phase annular region of annular centrifugal extractor. Ind. Eng. Chem. Res.
**2009**, 48, 37–46. [Google Scholar] [CrossRef] - Kolhe, N.S.; Mirage, Y.H.; Patwardhan, A.V.; Rathod, V.K.; Pandey, N.K.; Mudali, U.K.; Natarajan, R. CFD and experimental studies of single phase axial dispersion coefficient in pulsed sieve plate column. Chem. Eng. Res. Des.
**2011**, 89, 1909–1918. [Google Scholar] [CrossRef] - Fazli-Abukheyli, R.; Darvishi, P. Combination of axial dispersion and velocity profile in parallel tanks-in-series compartment model for prediction of residence time distribution in a wide range of non-ideal laminar flow regimes. Chem. Eng. Sci.
**2019**, 195, 531–540. [Google Scholar] [CrossRef] - Ding, S.; Liu, L.; Xu, J. A study of the determination of dimensionless number and its influence on the performance of a combination wastewater reactor. Procedia Environ. Sci.
**2013**, 18, 579–584. [Google Scholar] [CrossRef] - Elias, Y.; von Rohr, P.R. Axial dispersion, pressure drop and mass transfer comparison of small-scale structured reaction devices for hydrogenations. Chem. Eng. Process.
**2016**, 106, 1–12. [Google Scholar] [CrossRef] - Makinia, J.; Wells, S.A. Evaluation of empirical formulae for estimation of the longitudinal dispersion in activated sludge reactors. Water Res.
**2005**, 39, 1533–1542. [Google Scholar] [CrossRef] - Sharma, L.; Nigam, K.D.P.; Roy, S. Axial dispersion in single and multiphase flows in coiled geometries: Radioactive particle tracking experiments. Chem. Eng. Sci.
**2017**, 157, 116–126. [Google Scholar] [CrossRef] - Murphy, K.L.; Timpany, P.L. Design and analysis of mixing for an aeration tank. J. Sanit. Eng. Div.
**1967**, 93, 1–15. [Google Scholar] - Courant, R.; Hilbert, D. Methods of Mathematical Physics Volume II; Wiley-Interscience: Hoboken, NJ, USA, 1989. [Google Scholar]
- Thomas, H.A.; McKee, J.E. Longitudinal Mixing in Aeration Tank. Sew. Works J.
**1944**, 16, 42–55. [Google Scholar] - Borkowski, P. Adaptive system for steering a ship along the desired route. Mathematics
**2018**, 6, 196. [Google Scholar] [CrossRef] - Van der Laan, E.T. Notes on diffusion-type model for longitudinal mixing in flow. Chem. Eng. Sci.
**1958**, 7, 187. [Google Scholar]

t_{i} (min) | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 12 | 14 |

c_{i} (mg/L) | 0 | 1 | 5 | 8 | 10 | 8 | 6 | 4 | 3 | 2.2 | 1.5 | 0.6 | 0 |

i | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 11 | 12 | 13 |

${T}_{i}$ | 0 | 0.1979 | 0.3957 | 0.5936 | 0.7914 | 0.9893 | 1.1871 | 1.3850 | 1.5829 | 1.7807 | 1.9786 | 2.3743 | 2.77 |

i | ${\mathit{\mu}}_{\mathit{i}}$ | β_{i} | γ_{i} |
---|---|---|---|

1 | $2.114669$ | $0.162452$ | $1.468256$ |

2 | $4.525511$ | $0.215281$ | $1.848822$ |

3 | $7.239104$ | $0.195738$ | $1.350646$ |

4 | $10.133634$ | $0.163113$ | $0.445602$ |

$\mathit{i}$ | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
---|---|---|---|---|---|---|---|---|---|

${\left(\frac{D}{uL}\right)}_{i}$ | 0.13 | 0.14 | 0.15 | 0.16 | 0.17 | 0.18 | 0.19 | 0.2 | 0.21 |

$E{r}_{{\overline{c}}_{i}}\left({\left(\frac{D}{uL}\right)}_{i}\right)$ | 4.4090 | 3.9981 | 3.5411 | 3.2279 | 3.0203 | 2.8912 | 2.8213 | 2.7965 | 2.8061 |

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).