## Services on Demand

## Journal

## Article

## Indicators

- Cited by SciELO
- Access statistics

## Related links

- Cited by Google
- Similars in SciELO
- Similars in Google

## Share

## Brazilian Journal of Chemical Engineering

##
*Print version* ISSN 0104-6632

### Braz. J. Chem. Eng. vol.31 no.3 São Paulo July/Sept. 2014

#### https://doi.org/10.1590/0104-6632.20140313s00002625

**PROCESS SYSTEMS ENGINEERING **

**State estimation of chemical engineering systems tending to multiple solutions**

**N. P. G. Salau ^{I, *}; J. O. Trierweiler^{II}; A. R. Secchi^{III}**

^{I}Chemical Engineering Department, Universidade Federal de Santa Maria, UFSM, Av. Roraima 1000, Cidade Universitária, Bairro Camobi, CEP: 97105-900 Santa Maria, RS - Brazil. Phone: + (55) (55) 3220-8448, Fax: + (55) (55) 3220-8030 E-mail: ninasalau@ufsm.br

^{II}Chemical Engineering Department, Universidade Federal do Rio Grande do Sul, Porto Alegre - RS, Brasil

^{III}Chemical Engineering Program, COPPE, Universidade Federal do Rio de Janeiro, Rio de Janeiro - RJ, Brasil

**ABSTRACT**

A well-evaluated state covariance matrix avoids error propagation due to divergence issues and, thereby, it is crucial for a successful state estimator design. In this paper we investigate the performance of the state covariance matrices used in three unconstrained Extended Kalman Filter (EKF) formulations and one constrained EKF formulation (CEKF). As benchmark case studies we have chosen: a) a batch chemical reactor with reversible reactions whose system model and measurement are such that multiple states satisfy the equilibrium condition and b) a CSTR with exothermic irreversible reactions and cooling jacket energy balance whose nonlinear behavior includes multiple steady-states and limit cycles. The results have shown that CEKF is in general the best choice of EKF formulations (even if they are constrained with an ad hoc clipping strategy which avoids undesired states) for such case studies. Contrary to a clipped EKF formulation, CEKF incorporates constraints into an optimization problem, which minimizes the noise in a least square sense preventing a bad noise distribution. It is also shown that, although the Moving Horizon Estimation (MHE) provides greater robustness to a poor guess of the initial state, converging in less steps to the actual states, it is not justified for our examples due to the high additional computational effort.

**Keywords:** Nonlinear state estimation; State covariance matrix; Noise distribution; Multiple solutions.

**INTRODUCTION**

It is well known that a suitable design of state estimators requires a representative model for capturing the plant behavior and knowledge about the noise statistics, which are generally not known in practical applications (Valappil and Georgakis, 2000). However, some divergence issues such as numerical round-off, plant-model mismatch, and state unobservability also deserve special attention because they can lead to convergence failures (Brown and Hwang, 2012; Grewal and Andrews, 2008; Simon, 2006). Any state covariance matrix equation is composed of states, measurements, linear models and noise covariance statistics and, hence, all the mentioned divergence issues may increase the error propagation conveyed by this matrix.

Even with an accurate dynamic nonlinear process model, the estimation of noise covariance statistics is of utmost importance to obtain a well-computed state covariance matrix and thereby avoid the EKF failure. While the measurement noise covariance matrix can be directly derived from the measurement device, the determination of the process-noise covariance matrix is much less straightforward. The model accuracy is specified by tuning the Kalman filter that involves the selection of the process-noise covariance matrix. If this matrix is guessed low, the filter will believe the model excessively and will not use the on-line measurements properly to correct the states. This can lead to poor performance or even filter divergence. On the other hand, if the process-noise covariance matrix is guessed higher than the actual value, the state estimates will be noisy and uncertain, as this would lead to increased values of the state covariance matrix (Valappil and Georgakis, 2000; Salau *et al*., 2013). In a few words, choosing the right value of the tuning parameters is very important for successful application of EKF (Kalman, 1962).

In contrast to the measurement noise covariance, information about the uncertainty of the initial state guess cannot be measured experimentally. In general, the filter will not be fatally affected if the initial state guess is not close to the actual initial state, but convergence to the correct estimate may be slow. However, if the initial state estimate covariance matrix is chosen too small or optimistic while initial state guess and actual initial state differ considerably, the Kalman filter gain becomes small and the estimator relies on the model predictions more than it should.

Thus, subsequent measurements do not have the impact on the estimator that they need to have. The filter might then learn the wrong state too well and diverge (Jazwinski, 1970). The importance of choosing a consistent pair of an initial state estimate guess and initial state estimate covariance matrix is also emphasized by other authors (Valappil and Georgakis, 2000). For instance, Vachhani *et al*. (2004) emphasized the role of a consistent choice of initial state estimate covariance matrix and consequently the use of alternate initial state information for their simulations. Similarly, Prakash *et al*. (2010) silently adapt their initial state estimate to be more consistent with the original initial state estimate covariance matrix.

In this work, however, we will show that some EKF formulations can indeed converge, even in the presence of poor initial guesses, as well as a poor initial state estimate covariance matrix.

In the literature, several modified implementations of the EKF are presented in an effort to avoid the filter failure (Grewal and Andrews, 2008). The basic difference between these formulations is concerned with the state covariance matrix computation. We can also find in the literature some contributions which show examples and outline conditions where other methods perform better than the standard EKF, MHE for instance (Rao and Rawlings, 2002; Haseltine and Rawlings, 2005). Otherwise, these works address comparison issues such as accuracy and computational expense and no efforts are made to avoid divergence or even instability by EKF state covariance matrix computation and thereby improve the filter performance.

Due to the lack of literature concerning this problem, we investigate in this paper the performance of four EKF formulations, being three unconstrained: (i) the classical EKF formulation, called here discrete EKF (DEKF), (ii) EKF with the continuous Riccati equation (EKF-CRE) and (iii) reduced-rank extended Kalman filter (RREKF); and a constrained EKF (CEKF) formulation (Gesthuisen *et al*., 2001) to derive some results with two systems tending to multiple solutions giving insights into their numerical performance: a) a batch chemical reactor with reversible reactions whose system model and measurement are such that multiple states satisfy the equilibrium condition and b) a CSTR with exothermic irreversible reactions and cooling jacket energy balance whose nonlinear behavior includes multiple steady-states and limit cycles.

In the next section we briefly review the state estimation formulations of interest. Then, we present two systems tending to multiple solutions and outline conditions that lead a classical EKF formulation to converge to physically unrealizable equilibrium points and to an undesired steady-state. Then, we demonstrate the potentiality of EKF-CRE and RREKF to handle the divergence issues by state covariance matrix computation. Otherwise, both formulations may not prevent undesired states before eventually converging to the actual states. Although an ad hoc clipping strategy seems to be a reasonable solution to constrain the states, it disregards the assumption that the measurement noise is a Gaussian random noise.

Finally, we have shown that CEKF is in general the best alternative to EKF formulations due to the possibility of incorporating constraints into an optimization problem, hence preventing the estimator from converging to undesired states and from bad noise distribution. Furthermore, this technique demands a small computational effort and a performance comparable to the MHE.

**FORMULATION AND SOLUTION OF THE ESTIMATION PROBLEM**

Consider the following nonlinear dynamic and continuous-time system with discrete-time measurements to be used in the state estimation formulations

where *u* denotes the deterministic inputs, *x* denotes the states, and *y* denotes the measurements. The process-noise vector, , and the measurement-noise vector, , are assumed to be white Gaussian random processes with zero mean and covariance *Q* and *R*_{k}, respectively.

The system is linearized at each time step to obtain the local state-space matrices as below:

The filter algorithm is initialized by the initially expected state and state covariance:

**Unconstrained EKF Formulations**

The equations that compose the different steps in the considered unconstrained EKF formulations are given below.

State transition equation:

Kalman gain equation:

State update equation:

The basic difference of these formulations is concerned with the state covariance matrix computation, as described below.

**Discrete Extended Kalman Filter**

DEKF considers discrete-time dynamics and discrete-time measurements. This situation is often encountered in practice. Even if the underlying system dynamics are continuous-time, the state estimator usually needs to be implemented in a digital computer (Simon, 2006).

The nonlinear dynamic system in continuous time (cf. Eq. (1)) can be approximated by a linear dynamic system in continuous time expressed in a general form as a first-order vector of difference equation (Brown and Hwang, 2012; Grewal and Andrews, 2008; Simon, 2006).

where is the *state transition matrix* for the state at time *t*_{k} given as

The state covariance matrix transition equation of the DEKF is given by

and the state covariance matrix update equation can be represented by one of the following equations:

The first equation (Eq. (15)), called the Joseph form (Brown and Hwang, 2012; Grewal and Andrews, 2008; Simon, 2006), has somewhat better numerical behavior than the other forms in unusual numerical situations (Brown and Hwang, 2012) since it guarantees that will always be symmetric positive definite, as long as is symmetric positive definite, and *K*_{k} does not need to be optimal. The third equation for (Eq. (17)) is computationally simpler than the first equation, but its form does not guarantee symmetry or positive definiteness for (Simon, 2006). However, we have formulated the DEKF using the Eq. (17) for comparison with the reference of our first case study (Haseltine and Rawlings, 2005). Equation 16 is just the substitution of Eq. (10) into Eq. (17) and is rarely implemented as written above (Simon, 2006), but will be useful in our derivation of the information filter in the next sections.

**Reduced Rank Extended Kalman Filter**

The reduced rank extended Kalman filter (RREKF) expresses the transition state covariance matrix in terms of its dominant eigenvectors. Neglecting non-dominant eigenvectors of the covariance matrix implies that the dimensionality of the confidence region is reduced. Measurement updates therefore have no effect on the directions of non-dominant eigenvectors. If the actual state does not leave the attractor, the RREKF is advantageous since measurement updates are omitted in those directions of strongest attraction. Furthermore, this EKF modification is more robust against improper initialization, as claimed by Pham *et al* (1998).

Although the RREKF is not applied or at least barely applied to chemical processes, we have selected this formulation to show the benefits of a lower rank covariance approximation of the transition state covariance matrix (Eq. (14)) in the examples we are interested in. Except for this modification, RREKF presents the same formulation as for DEKF.

We briefly outline the implementation and implications of the low-rank covariance approximation. Let comprise the eigenvectors of as columns and the corresponding eigenvalues on the diagonal. As the transition state covariance matrix (Eq. (14)) is symmetric, it can be decomposed as

For , the rank approximation of the covariance matrix is given by

Geometrically, the lower-rank approximation is the orthogonal projection of the covariance ellipse or (hyper-)ellipsoid onto its *q* most dominant axes.

**Extended Kalman Filter with the Continuous-Time Riccati Equation**

In this section we introduce an EKF formulation referred to as EKF with the continuous-time Riccati equation (EKF-CRE).

Before we define the EKF-CRE, a brief review about an alternate form for the EKF which applies the discrete-time Riccati equation (DRE) is necessary (Simon, 2006). This approximation referred to here as EKF-DRE is computationally attractive because the state covariance matrix, *P*, is evaluated just once at each time step, i. e. the state covariance matrix transition and update steps are carried out simultaneously.

Considering Eq. (16) for the state covariance matrix update, rewritten here for convenience:

Applying Eq. (14) considering the next time steps, we can write the state covariance matrix transition as follows:

Substituting Eq. (20) into Eq. (21) gives

Rearranging Eq. (22) and assuming that , because the state covariance matrix transition and update steps are carried out simultaneously, we obtain the equation for the one-step state covariance matrix equation

Now, rewriting Eq. (23), considering one time step before, we arrive to the formulation called DRE

Note that the Kalman gain (Eq. (10)) requires and this term is not computed in this alternate EKF formulation. Thus, Eq. (10) can be modified by assuming that , resulting in

In spite of applying the one-step state covariance matrix equation and using a different equation for the Kalman gain, the EKF-DRE and the DEKF result in identical results and estimation-error covariances (Simon, 2006). To avoid a redundancy, the EKF-DRE will be not compared in our examples.

Alternatively, following the same steps, we can write a similar expression for the continuous case, which can be written as follows:

The equation above is known as the continuous-time Ricatti equation (CRE). It is equivalent to the state covariance matrix transition equation of the continuous-time EKF found in the literature (Brown and Hwang, 2012; Grewal and Andrews, 2008; Simon, 2006). The basic difference between both formulations is concerned with the Kalman gain equation. For the continuous-time EKF, the Kalman gain equation considers that the measurement values remain constant during the entire time interval, which is suitable just for small sampling times. Furthermore, the EKF-CRE considers continuous-time dynamics and discrete observation for update, i. e. discrete-time measurements.

**Moving Horizon Estimator (MHE)**

Before explaining the CEKF formulation, the basic aspects of the MHE (Muske and Rawlings, 1994; Robertson *et al*., 1996; Rao *et al*., 2003) are presented. As paper of Haseltine and Rawlings (2005), we examine the performance of MHE with local optimization and an arrival cost approximated with a smoothing update (Tenny and Rawlings, 2002). According to Tenny and Rawlings (2002) the smoothing scheme is superior to the filtering scheme because the filtering scheme induces oscillations in the state estimates through the unnecessary propagation of initial error. Further details regarding this MHE scheme can be found in the mentioned manuscript.

The basic idea of MHE is to proceed with state estimation by using only the most recent *N+*1 measurements, where *N* is the time horizon size.

The moving horizon approximation of the objective function is given by

subject to the equality constraints

and the inequality constraints

From optimization, measurements and states are updated as follows

Rao *et al*. (2003) suggested computing the state covariance matrix equation (Eq. (31)) recursively using the discrete Riccati equation. One obtains this result by deriving the deterministic Kalman filter using forward dynamic programming (Cox, 1964).

Note that the equation above is the same as for Eq. (24) applied to each horizon step.

We solve Eq. (27) using sequential quadratic programming (SQP) as implemented in the medium-scale algorithm of the *fmincon* function in MatLab. For the successive integration of Eq. (30) we use DASSLC (Differential-Algebraic System Solver in C) which does the multirate integration of systems of differential-algebraic equations (DAEs). The integration algorithm used in DASSLC is an extension of the DASSL code of Petzold (1983) to address high-index and large-scale problems, and the setup algorithm used in DASSLC is based on the DAWRS code (Secchi *et al*., 1991), a package to solve DAEs on parallel machines.

**Constrained EKF**

CEKF follows from the MHE when the horizon length is set to zero (Gesthuisen *et al*., 2001). Zero length implies that ODEs are not considered in the optimization problem, which simplifies the complexity of solving a nonlinear dynamic optimization problem.

Setting (*N=*0) in the MHE optimization problem (Eq. (27)), the resulting formulation is exactly the CEKF formulation problem:

subject to the equality constraints:

and inequality constraints:

From optimization, measurements and states are updated as follows:

If the measurement equation is linear, the resulting problem is a quadratic program (QP), which can be solved with small computational effort (Gesthuisen *et al*., 2001). We use the MatLab's *quadprog* function to solve the quadratic programming problem as given below.

where:

subject to:

In contrast to the unconstrained EKF formulations, state estimation with CEKF is formulated as an optimization problem, so that constraints on state variables can be incorporated in this problem. Furthermore, it is not necessary to compute the Kalman gain (Eq. (10)) and therefore the state covariance matrix *P* is computed using the discrete Riccati equation as in Eq. (24).

Although the CEKF belongs to a family of MHE estimators, we choose to call this formulation here CEKF, as introduced by Gesthuisen *et al*. (2001). Note that CEKF is similar to EKF because the latter also takes into account a zero-length horizon in the updating stage, i. e. the current state is estimated only using the current measurement on the horizon and it is equivalent to zero horizon length (*N=*0). In the absence of any constraints and for low process-noise uncertainties, these formulations are similar.

Afterwards, Vachhani *et al*. (2005) have proposed essentially the same formulation as for CEKF with a different name: Recursive Nonlinear Dynamic Data Reconciliation (RNDDR). However, we have decided to adopt the earlier nomenclature.

**EXAMPLES OF DIVERGENCE ISSUES BY STATE COVARIANCE MATRIX COMPUTATION**

In this section, we outline the conditions that generate divergence issues by state covariance matrix computation in two chemical engineering systems. The first one is a batch reactor with reversible reactions whose system model and measurement are such that multiple states satisfy the equilibrium condition. The results obtained with this system (Figures 1 to 6) have been published in our previous work (Salau *et al*., 2009).

The second one is a CSTR with exothermic irreversible reactions and cooling jacket energy balance whose nonlinear behavior includes multiple steady-states and limit cycles.

**System Tending to Multiple Equilibrium Points**

The first example was introduced by Haseltine and Rawlings (2005) and consists of a gas-phase batch reactor with multiple equilibrium points considering the following reversible reactions:

The model equations and the model parameters are given below.

with

Multiple states satisfy the equilibrium condition for a given measurement, which in this case is the system pressure at the equilibrium, evaluated by the following equation:

Table 1 shows the possible theoretical solutions, without measurement or state noise, at the equilibrium pressure given by the initial state:

Note that only the solution *EQ*1 has physically realizable states (non-negative concentrations).

The state estimator parameters and the poor initial guess* *used for this example were obtained from Haseltine and Rawlings (2005):

According to the authors, EKF fails in this example because the system model and measurement are such that multiple states satisfy the equilibrium condition and is given a poor initial guess of the initial state for the estimator. Nonetheless, we cannot assert that EKF fails because two equilibrium points (*EQ*1 and *EQ*2 in Table 1) are possible due to the poor and incoherent guesses of the initial state and its covariance matrix used in this example, as discussed in the next section.

**RESULTS AND ANALYSES**

**Measurement Noise Perturbation**

First, uniform and normally distributed random measurement noises were simulated. Either solution *EQ*1 or solution *EQ*2 (i. e. the inconsistent solution) is obtained in accordance with the set of random measurement noise employed in DEKF. Thus, we have chosen a set of random measurement noise which leads DEKF to converge to the solution *EQ*1 and added a noise measurement perturbation in order to lead it to converge to solution *EQ*2. As can be seen in Figure 1, a noise measurement perturbation of 0.754 *atm* at *t*=0.5 *min* changes the estimated states trajectory from solution *EQ*1 to solution *EQ*2*.*

We have also found that other divergence issues can lead DEKF to converge to solution *EQ*2, such as low-precision arithmetic, numeric system linearization by finite differences, a poor guess of *P*_{0} and incorrect values of the tuning parameters: *Q* and *R*.

**Comparison Between Unconstrained EKF Formulations**

RREKF and EKF-CRE were applied to the batch reactor in order to prevent the state estimator from converging to *EQ*2. In spite of a measurement noise perturbation, both unconstrained EKF formulations converged to solution *EQ*1, as shown in Figure 2.

RREKF disregards non-dominant eigenvectors, which implies zero variance in the respective directions, and no effect of measurement updates. For this example, we use a rank-two approximation, so the confidence region is an ellipse instead of a (hyper-) ellipsoid, reducing the dimensionality of the confidence region. After initial oscillations, the estimated states converge to the actual states. At the equilibrium point , we observe that the eigenvector corresponding to the smallest eigenvalue or the non-dominant eigenvector is orthogonal to the tangent of the equilibrium curve (the scalar product approaches zero). This eigenvector is not considered in the rank-two covariance approximation and thus the filter applies no correction in this direction orthogonal to the attractor.

EKF-CRE is not subject to errors due to model discretization. The state covariance matrix computed by EKF-CRE in a single stage (Eq. (26)) guarantees that it will always be symmetric positive definite, since it is derived from the called Joseph form (Brown and Hwang, 2012; Grewal and Andrews, 2008; Simon, 2006), as discussed previously in the present article. Besides, the state covariance matrix of EKF-CRE presents a smaller condition number, i. e. it is less sensitive to perturbations than the states covariance matrices computed by DEKF (eqs 14 and 17), as shown in Figure 3. The mentioned advantages of EKF-CRE over DEKF justify the convergence of this formulation to *EQ*1, even with a measurement noise perturbation.

Although EKF-CRE and RREKF prevented physically unrealizable states at the final batch time, physically unrealizable states (negative concentrations) were unavoidable during the batch. This fact is justified by the poor guesses of the initial state and its state covariance matrix.

We also considered that the unconstrained formulations fail when they converge to negative concentrations (*EQ2*) because they do not converge to the given plant measurements (that are consistent with positive concentrations - *EQ1*), but instead, converge to physically unrealizable states, a fact that could not be considered as a success from the chemical engineering point of view. In the next section it will be shown that constrained formulations can avoid such problems.

**Comparison Between Clipped DEKF and CEKF**

To prevent physically unrealizable states, we have constrained DEKF with an ad hoc clipping strategy (Haseltine and Rawlings, 2005) in which negative update values of the state are set to zero (i. e. if ).

The comparison between clipped DEKF and CEKF performances is shown in Figure 4. Note that before eventually converging to the actual states, the pressure filtered by the clipped DEKF is quite larger than the measured one. However, the clipped DEKF in our work presented a better performance when compared with the one presented by Haseltine and Rawlings (2005). In their paper, the clipped DEKF drives the predicted pressure 3 orders of magnitude larger than the measured pressure before eventually converging to the actual states on a longer time scale (1 order larger than the convergence time obtained in our work). This occurred because the authors overestimated *P*_{0} (cf. Eq. (47)). While the initial guess for the concentration of specie C is too far from its actual value, the same is not valid for the concentration of specie B. Since the authors have chosen a *P*_{0} with the same weight for all initial guesses (diagonal elements of *P*_{0}), this weight should be balanced between the initial guesses. We have selected, hence, a lower initial guess of the state covariance matrix:

Thus, due to a lower initial guess of the state covariance matrix, the clipped DEKF is shown as a reasonable alternative to prevent physically unrealizable states.

A clipped DEKF, however, disregards the assumption that is a Gaussian random noise and does not let the Kalman gain properly distributes the measurement residuals throughout the estimated state and, thereby, corrects them. Another disadvantage of such an approach is that clipped state estimates violate the process model, as the negative state estimates were obtained from the model solution.

On the other side, CEKF swiftly converges to the actual states and minimizes and in a least square sense, incorporating constraints into an optimization problem, which prevents bad noise distribution. Figure 5 illustrates the distribution of for the clipped DEKF and the CEKF. It means that the main difference between the clipped DEKF and CEKF is that the CEKF preserves the gaussianity, which is one of the main assumptions of the Kalman Filter Approach.

It is well known that the quality of the MHE estimates is a function of the estimation horizon (Haseltine and Rawlings, 1970). Thus, we have enlarged the MHE estimation window to *N*=2 and compared its performance with the CEKF, as can be seen in Figure 6. In this comparison, CEKF and MHE are computed recursively using the discrete Riccati equation (DRE) (Eq. (31)).

The Integral Time Absolute Error (ITAE) has been chosen as the performance criterion to compare CEKF and MHE. The results obtained with ITAE, as well as the computational expenses of both filters, are shown in Table 2.

Even with the MHE considering a size horizon of 4, the average CPU time per time step is smaller than the sampling period (0.25 min). Thus, both CEKF and MHE are suited for this case study.

However, MHE presented a computational effort 2 orders of magnitude larger than CEKF and, hence, it is not justified for this process as the improvement of estimation results is marginal: from 12 to 20% considering a size horizon of 2 and from 15 to 20% considering a size horizon of 4. Therefore, due to the small computational effort of CEKF, the possibility of incorporating constraints into an optimization problem, and the estimation results comparable to MHE, this approach seems to be the best choice for such a case study.

**System Tending to Multiple Steady-States and Limit Cycles**

As a benchmark example, we have chosen a CSTR, as introduced by Torres and Tlacuahuac (2000). The following two exothermic irreversible first-order reactions in series take place in the reactor:

The reactor volume and physical parameters are assumed to remain constant; perfect mixing is also assumed. In addition the dynamics of the cooling jacket are taken into account.

The dimensionless model equations are given below and the dimensionless parameters are defined in Table 3. More details on the model can be found in Torres and Tlacuahuac (2000).

where *x*_{1} is the dimensionless concentration of reactant *A*, *x*_{2} is the dimensionless concentration of reactant *B*, *x*_{3} is the dimensionless reactor temperature, and *x*_{4} is the dimensionless cooling jacket temperature. The dimensionless cooling water volumetric flowrate, *q*_{c}, is the manipulated variable.

Torres and Tlacuahuac (2000) have analyzed input/output multiplicities of the full model using *q*_{c} as continuation parameter. In the bifurcation diagram of Figure 7, five steady-states and a bifurcation point were observed for *x*_{3} when *q*_{c}=2.3.

The CSTR model parameter values used for generating Figure 7 are shown in Table 4.

We define the state and measurements to be:

We consider state estimators with the following parameters:

**RESULTS AND ANALYSES**

In this section we analyze the state estimator performances considering the operating region in which five steady-states and a Hopf bifurcation point were observed for *x*_{3} when *q*_{c}=2.3 (cf. Figure 7).

**Comparison Between Unconstrained EKF Formulations**

The poor initial guess given by Eq. (64) was selected because it leads DEKF to converge to an undesirable steady-state. In the absence of any measurement noise perturbation, RREKF and EKF-CRE were also applied to the CSTR. As shown in Figure 8, both unconstrained formulations converged to the actual steady-state, which in fact is a limit cycle. Likewise for the first example, RREKF presents a slower convergence to the actual steady-state when compared to EKF-CRE.

As mentioned before, RREKF disregards non-dominant eigenvectors. For this example, we use a rank-one approximation and three non-dominant eigenvectors were eliminated:

At the equilibrium point , we observe that the eigenvectors and are bi-orthogonal to the tangent of the equilibrium curve (the scalar product is negative). These three eigenvectors are not considered in the rank-one covariance approximation and thus the filter applies no correction in their directions.

Therefore, the RREKF applies correction only in the dominant eigenvector direction: . However, this estimator converges slowly to the actual steady-state for this example.

Again, the state covariance matrix computed by EKF-CRE in a single stage presents a smaller condition number than the state covariance matrices computed by DEKF, as shown in Figure 9.

**Comparison Between Clipped DEKF and CEKF**

We have constrained DEKF with an ad hoc clipping strategy in order to prevent an undesirable steady-state. Nonetheless, the strategy of resetting to zero the negative update values did not perform well for this example and the clipped DEKF converged to the same undesirable steady-state as for the DEKF (cf. Figure 8). Thus, we selected stricter constraints for update values of the estimated states *x*_{3} and *x*_{4}:

if

if

The strategy of resetting to zero the negative update values was maintained for the measured states *x*_{1} and *x*_{2}*.*

The comparison between clipped DEKF and CEKF performances are shown in Figure 10.

Although both formulations swiftly converged to the actual states, the CEKF performed slightly better than the clipped DEKF for the CSTR case, as shown in Figure 10.

Finally, we also enlarged the MHE estimation window to *N*=2 for the CSTR case. In addition, the state covariance matrix of MHE was computed recursively using either the discrete Riccati equation (MHE-DRE) (Eq. (31) or the continuous Riccati equation (MHE-CRE) (Eq. (26) applied to each horizon step). The comparison between the CEKF, MHE with DRE (*N*=2) and MHE (*N*=2) with CRE performances is shown in Figure 11.

According to Figure 11, MHE is superior to CEKF to swiftly converge to the actual states. Since CRE is not subject to errors due to model discretization, it is not surprising that MHE-CRE presents more accurate estimates of states.

The comparisons between the Integral Time Absolute Error (ITAE) and the computational expenses of CEKF, MHE-DRE and MHE-CRE are presented in Table 5. Such values were calculated for 20 time units, while Figures 10 and 11 were illustrated up to 2 time units for better visualization of the convergence. In this case, the improvement of the estimation when MHE is implemented is expressive: around 80% considering a size horizon of 2 in MHE-DRE and around 90% considering a size horizon of 2 in MHE-CRE. However, the improvement of estimation when MHE is implemented is reflected in a higher computational effort: 3 orders of magnitude larger.

In the case study reference (Torres and Tlacuahuac, 2000), the variables are dimensionless and nothing is said about the original units. However, in order to estimate the non-measured states, we suppose a sampling time of 0.1. If this value were in hours, which is proper to obtain concentrations measurements, both CEKF and MHE are suited for this case study.

It has been seen that MHE provides improved state estimation and presents greater robustness to a poor guess of the initial state. However, after converging to the actual states, both MHE and CEKF perform equally accurately and, therefore, the use of MHE becomes needless for the systems studied in this work.

In practice, however, chemical engineering systems are frequently subject to unexpected process disturbances. Because MHE employs a trajectory of measurements during the state estimation, it shall present a better performance when compared to CEKF to handle with such problem.

**CONCLUSION**

This paper outlines the performance of the state covariance matrices used in three unconstrained Extended Kalman Filter (EKF) formulations and one constrained EKF formulation (CEKF), for two chemical engineering examples tending to multiple solutions.

The first example is a batch reactor with reversible reactions whose system model and measurement are such that multiple states satisfy the equilibrium condition. With a measurement noise perturbation, we outline a condition that led a classical EKF formulation (DEKF) to converge to physically unrealizable or undesired states. According to our results, EKF-CRE and RREKF are more numerically robust in computing the state covariance matrix than DEKF. As both formulations avoided an increase in error propagation due to a measurement noise perturbation, they were able to converge to the actual steady-state. Thus, a suitable choice of the EKF formulation based on the state covariance matrix equation is essential to prevent physically unrealizable states. However, due to the poor guesses of the initial state and its state covariance matrix for the batch reactor, EKF-CRE and RREKF did not prevent physically unrealizable states during the batch. The clipped DEKF has been seen as a reasonable alternative to prevent physically unrealizable states, although this approach converges slowly to the actual states and disregards the assumption that the measurement noise is a Gaussian random noise and that clipped state estimates violate the process model. On the other hand, the CEKF can be seen as the best technique for systems with such behavior due to the possibility of incorporating constraints into an optimization problem minimizing the noise in a least square sense, preventing bad noise distribution. Besides, it is simpler, computationally less demanding than the MHE, and has comparable performance.

The second example is a CSTR with exothermic irreversible reactions and cooling jacket, whose nonlinear behavior includes multiple steady-states and limit cycles. The results for the CSTR demonstrate that, similar to the batch reactor case, EKF-CRE and RREKF converged to the actual steady-state. However, the clipped DEKF and CEKF presented faster convergence to the actual states, and the CEKF performed slightly better than the clipped DEKF for this example. Contrary to the batch reactor, MHE presented a superior performance for the CSTR case. Further, we demonstrated that MHE performs better when the state covariance matrix is computed recursively using the continuous-time Riccati equation.

Since MHE employs a trajectory of measurements as opposed to measurements at only a single time, it is better suited than the CEKF (MHE with a zero horizon length) to handle a poor guess of the initial state and unexpected process disturbances. In these circumstances, it can be possible for the CEKF to fail to converge swiftly to the actual states. Hence, an adaptive strategy to select the MHE horizon length shall be considered for future work. For this purpose, the state estimation quality is evaluated through sensitivity analyses to trigger the horizon length augment when necessary. This approach aims to handle effectively the divergence issue while requiring only a minimum size of horizon length, reducing the online computational expense.

**ACKNOWLEDGMENT**

We gratefully acknowledge Prof. James Rawlings and Prof. Fernando Lima for the fruitful discussions. The first author also would like to thank the financial support from the German Academic Exchange Service (DAAD) and from PETROBRAS S.A.

**REFERENCES**

Brown, R. G. and Hwang, P. Y. C., Introduction to Random Signals and Applied Kalman Filtering: With MATLAB Exercises and Solutions, 4th Edition. John Wiley & Sons, New York (2012). [ Links ]

Cox, H., On the estimation of state variables and parameters for noisy dynamic systems. IEEE Transactions on Automatic Control, 9(1), 5-12 (1964). [ Links ]

Gesthuisen, R., Klatt, K. -U., Engell, S., Optimization-based State Estimation: A comparative study for the batch polycondensation. European Control Conference (ECC) (2001). [ Links ]

Grewal, M. S. and Andrews, A. P., Kalman Filtering: Theory and Practice Using MatLab, 3rd Edition, John Wiley & Sons, New Jersey (2008). [ Links ]

Haseltine, E. L. and Rawlings, J. B., Critical evaluation of extended Kalman filtering and moving-horizon estimation. Ind. Eng. Chem. Res., 44, 2451-2460 (2005). [ Links ]

Jazwinski, A. H., Stochastic Processes and Filtering Theory. Academic Press, New York (1970). [ Links ]

Kalman, R., Canonical structure of linear dynamical systems. Proc. Natl. Acad. Sci., 48, 596-600 (1962). [ Links ]

Muske, K. R. and Rawlings, J. B., Nonlinear moving horizon state estimation, methods of model based process control. NATO ASI Series, 293, 349-365 (1994). [ Links ]

Petzold, L. R., DASSL: Differential Algebraic System Solver. Technical Report of Sandia National Laboratories. Livermore (1983). [ Links ]

Pham, D. T., Verron, J. and Roubaud, M. C., A Singular evolutive extended Kalman filter for data assimilation in oceanography. Journal of Marine Systems, 16(3-4), 323-340 (1998). [ Links ]

Prakash, J., Patwardhan, S. C. and Shah, S. L., Constrained nonlinear state estimation using ensemble Kalman filters. Ind. Eng. Chem. Res., 49, 2242-2253 (2010). [ Links ]

Rao, C. V., Rawlings, J. B., Constrained process monitoring - moving horizon approach. AIChE Journal, 48(1), 97-109 (2002). [ Links ]

Rao, C. V., Rawlings, J. B. and Mayne, D. Q., Constrained state estimation for nonlinear discrete-time systems: Stability and moving horizon approximations. IEEE Transactions on Automatic Control, 48(2), 246-258, (2003). [ Links ]

Robertson, D. G., Lee, J. H. and Rawlings, J. B., A moving horizon-based approach for least-squares estimation. AIChE Journal, 42(8), 2209-2224 (1996). [ Links ]

Salau, N. P. G., Trierweiler, J. O. and Secchi, A. R., Numerical pitfalls by state covariance computation. Computer Aided Chemical Engineering, 27, 1215-1220 (2009). [ Links ]

Salau, N. P. G., Trierweiler, J. O. and Secchi, A. R., Observability analysis and model formulation for nonlinear state estimation. Accepted for Publication in Applied Mathematical Modeling (2013). [ Links ]

Secchi, A. R., Morari, M., Biscaia Jr., E. C., DAWRS: A Differential - Algebraic System Solver by the Waveform Relaxation Method. The Sixth Distributed Memory Computing Conference (DMCC6), Portland (1991). [ Links ]

Simon, D., Optimal State Estimation: Kalman, H-infinity, and Nonlinear Approaches. John Wiley & Sons, New Jersey (2006). [ Links ]

Tenny, M. J. and Rawlings, J. B., Efficient moving horizon estimation and nonlinear model predictive control. American Control Conference, Anchorage (2002). [ Links ]

Torres, A. E. G. and Tlacuahuac, A. F., Effect of process modeling on the nonlinear behaviour of a CSTR reactions AàBàC. Chemical Engineering Journal, 77, 153-164 (2000). [ Links ]

Vachhani, P., Narasimhan, S. and Rengaswamy, R., Recursive state estimation in nonlinear processes. American Control Conference, Boston (2004). [ Links ]

Vachhani, P., Rengaswamy, R., Gangwal, V. and Narasimhan S., Recursive estimation in constrained nonlinear dynamical systems. AIChE Journal, 51(3), 946-959 (2005). [ Links ]

Valappil, J. and Georgakis, C., Systematic estimation of state noise statistics for extended Kalman filters. AIChE Journal, 46(2), 292-398 (2000). [ Links ]

Submitted: March 26, 2013

Revised: September 5, 2013

Accepted: September 27, 2013

* To whom correspondence should be addressed