Journal of Control Science and Engineering

Volume 2008 (2008), Article ID 154956, 11 pages

http://dx.doi.org/10.1155/2008/154956

## Empirical Reduced-Order Modeling for Boundary Feedback Flow Control

^{1}Department of Electrical Engineering and Computer Science, University of Tennessee, 1508 Middle Drive, Knoxville, TN 37996, USA^{2}Performance Assessment and Decision Analysis Department, Carlsbad Program Group, Sandia National Laboratories, 4100 National Parks Highway MS 1395, Carlsbad, NM 88220, USA^{3}Air Force Research Laboratory, Wright-Patterson Air Force Base, Dayton, OH 45433, USA

Received 27 December 2007; Revised 4 November 2008; Accepted 24 December 2008

Academic Editor: Onur Toker

Copyright © 2008 Seddik M. Djouadi et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

This paper deals with the practical and theoretical implications of model reduction for aerodynamic flow-based control problems. Various aspects of model reduction are discussed that apply to partial differential equation- (PDE-) based models in general. Specifically, the proper orthogonal decomposition (POD) of a high dimension system as well as frequency domain identification methods are discussed for initial model construction. Projections on the POD basis give a nonlinear Galerkin model. Then, a model reduction method based on empirical balanced truncation is developed and applied to the Galerkin model. The rationale for doing so is that linear subspace approximations to exact submanifolds associated with nonlinear controllability and observability require only standard matrix manipulations utilizing simulation/experimental data. The proposed method uses a chirp signal as input to produce the output in the eigensystem realization algorithm (ERA). This method estimates the system's Markov parameters that accurately reproduce the output. Balanced truncation is used to show that model reduction is still effective on ERA produced approximated systems. The method is applied to a prototype convective flow on obstacle geometry. An feedback flow controller is designed based on the reduced model to achieve tracking and then applied to the full-order model with excellent performance.

#### 1. Introduction

Recently there has been significant interest in model reduction for the purpose of control design [1–9]. One such application of reduced-order modeling is control design in the context of aerodynamic flow. Aerodynamic flow control is a research area of great interest to the air force and the fluid mechanics community. Presently considerable research efforts are working with feedback control law design for systems described by PDEs that need a very large number of states to accurately simulate their characteristics. However, recent advances in the design of actuators and sensors can be leveraged for better system control only if the control design methods provide a reliable low-order controller [10]. Additionally, simulation, and experimental diagnostics are making applications such as the suppression of acoustic tones in cavities, separation control for high lift, and trajectory control without the need to move hinged surfaces a possibility [11]. However, these applications require the integration of feedback control because of the need for robustness to flight condition and vehicle attitude, precision tracking, overcoming low-fidelity models, or moving a system away from a stable solution or limit cycle as efficiently as possible [11].

The traditional systematic development of feedback control laws for these systems subject to a large number of states is currently an intractable problem. Feedback control strategies offer the possibility to improve performance and reduce control power through the control of unstable structures in the flow field. Reduced models are important for the design of feedback control laws, which rely on models that capture the relevant dynamics of the input-output system and are amenable to control design.

Unfortunately, it is very difficult to create models that capture the relevant dynamics of the input-output system. For example, computational fluid dynamics simulations can provide good solutions to a discretized version of the Navier-Stokes equation [1, 12]. However, accurate simulations for simple shapes such as two-dimensional airfoils, or complex shapes, such as a full vehicle, require several thousands to millions of states. Therefore, the simulation results are not directly useful for control design [11]. Complexity in the model is a legitimate need. The large number of states is necessary to capture important flow features that occur at extremely small spatial scales. Although these small flow features might seem insignificant, if they are not captured, it is not possible to analyze if they are necessary in securing the closed-loop system’s overall stability [10].

POD has been extensively investigated in distributed parameters systems due to its order reduction capability [7–9, 13, 14], and balanced truncation, which is a simple yet efficient model reduction technique widely used in reducing model orders of high-order linear systems [15–17]. POD models of only a few dozen states can often accurately capture the input-output behavior of systems that have full-order system models of thousands of states [11]. In addition to using the POD method in conjunction with model reduction techniques, the idea of using empirical gramians is growing in popularity for use in an approximate balanced truncation [4, 18–20]. Further, some work has been done on finding nonlinear empirical Gramians for balanced truncation [20, 21].

In fluid flow configurations it is not uncommon for discretized flow models to describe thousands to millions of state variables, for example, if one uses a linear quadratic regulator (LQR) control formulation, roughly Riccati unknowns need to be calculated for a discretized flow model describing states. Existing computing power and computational algorithms are not capable of solving an LQR problem of such large dimension. For dynamical models that are very large scale, such as those describing fluid flow configurations, it is apparent that the order of the system must be reduced prior to control law design [22]. This prevents us from using closed-loop model reduction strategy wherein the system is part of a closed-loop system with a controller [23].

In the area of fluid mechanics controls must often be fixed to the boundary of the problem geometry. For example, control of flow separation over an airfoil requires that actuation and sensing be done on the airfoil surface [11]. The problem geometry used for this project is one example of a case where control is restricted to the boundaries by physical necessity.

The paper is organized as follows. In Section 2, we introduce a prototype flow problem geometry that is used to apply the proposed order reduction techniques. Section 3 introduces the eigensystem realization algorithm (ERA) to identify the Markov parameters of the system, and as a product of the empirical gramians. Section 4 introduces empirical balanced truncation. This method is based on approximate (empirical) controllability and observability gramians and uses only a single simulation/experimental test. In Section 5, empirical balanced truncation and the ERA algorithm are applied to the Galerkin model, and numerical results are provided to show the effectiveness of the proposed method. In Section 6, an controller based on the empirical reduced model and which achieves tracking is discussed. The responses of the controlled closed loop on the full-order model are presented and show that the controller achieves good tracking performances despite being designed on a much lower model than the original model. Section 7 contains concluding remarks.

#### 2. Problem Geometry

The specific problem geometry considered is shown in Figure 1. The idea and methods presented here could be modified to apply to a different geometry or obstacle shape. The problem statement with its corresponding boundary conditions and governing equations was taken from [10]. A realistic example of this geometry in an aerodynamic application would be a payload hatch open during flight with actuator control only on the boundary. Let be the region defined by Let be the region defined by Then the problem domain is given by In this problem setup, is an obstacle.

The system dynamics that act within the problem domain are described by the two-dimensional (2D) Burgers’ equation [10] where the form of is In this case, the value for is equal to 1 and is equal to 0. The value is similar to the Reynolds number used in the Navier-Stokes equation. This parameter controls how much nonlinearity is present in the problem. The value used is 300, a small “Reynolds number,” but it still allows for the nonlinearity to show in the problem.

Dirichlet boundary conditions located on the obstacle top and bottom are denoted by A Dirichlet boundary condition is a first-type boundary condition that specifies the values of the solution defined by on a domain boundary [24]. The form of the boundary condition is The boundary conditions on the top and bottom are described by the following: Here and are the control inputs on the top and bottom boundaries, respectively; the spatial functions and describe the spatial effect that the controls have on the top and bottom boundaries.

The boundary condition on the airflow
intake side is and it is parabolic in nature. The
airflow outtake side has a Neumann boundary condition that has the form [24] On all of the remaining boundaries of is set equal to 0 for all values of *t.* Finally, the initial conditions for
the interior are given by A
numerical solution was found by simulation using a uniformly spaced grid. The
resulting system model contains a little more than 2000 states [11]. The POD
model construction is based on the total energy captured. A condition that 99.9% of system energy must
be captured was used for determining how many system modes were retained. This condition was met by a 40 POD basis.
Although this is a major reduction from the numerical solution, it will be
shown that important system dynamics can be retained with even lower state
number system models.

The general approach of POD is to construct a series of solution “snapshots.” These snapshots are generated by numerical simulations of the governing system equation(s) with a variety of input equations [6, 25]. The snapshots are needed to generate the correlation matrix that is used to find the POD basis. It is important to choose relevant input signals for the numerically simulated system. Further, these inputs should be similar to the expected inputs of the real system. The inputs used for system identification are of the form [10] where the values for are and and the range for is 0 to 10 seconds with a sample every 50 milliseconds. The squelch signal for all three values of is shown in Figure 2.

The numerical simulation was
performed to create the ensemble of solution snapshots [10]. The value for *M* (the number of snapshots) must be greater than the number of
modes that one will choose for the approximated system model. For a good
representation *M* should be much
larger than the desired size for the POD basis [10].

The solution to the PDE is assumed be finite energy, that is, belongs to the Hilbert space The solution can be approximated as where the 's are time varying coefficients that multiply the spatial functions 's and where and are the standard Hilbert spaces of absolutely square integrable functions defined, respectively, on the time interval and spatial domain The approximation (9) can be as accurate as desired since the tensor space is dense in [18].

Any basis for can be used to construct the approximation of the solution Here we use the POD basis since it is optimal in the following sense [2]: In [2], it is shown that solving the above optimization reduces to the usual POD optimality in the average kinetic energy sense (15) discussed below [25].

The POD basis function is chosen to maximize the average projection
of the member onto suitably normalized [25] To construct numerically the POD
basis
we build the correlation matrix *L*, of
size composed of the inner products
of the snapshots. In we have where denotes complex conjugate transpose.

A spectral
decomposition of the matrix *L* is
performed. The *n* largest eigenvalues of the matrix *L* are found and placed in descending order. Then the set of
eigenvectors are identified to be

The resulting orthonormal POD basis of dimension *n* can be constructed using the information found from the
correlation matrix *L*. First, the
eigenvectors of *L* are weighted by
their corresponding eigenvalues and normalized according to [10] Then, the POD basis set is formed according to with being the *j*th component of the eigenvector
Solving (17) gives 's,
which constitute the POD basis of dimension *n*.

The governing equation is projected onto the POD basis. The projection is accomplished via a Galerkin type projection and results in a system of ordinary differential equations (ODEs). The Galerkin projection results in only a weak solution to the PDE. However, this weak solution with finite difference approximations of the boundary conditions eventually leads to a nonlinear temporal model for the temporal or POD coefficients [10].

Projecting (1) onto the POD basis yields
[10] where the first term on the right-hand
side is The Neumann boundary condition forces the
portion of the boundary integral over to along to be 0, that is, The second boundary integral is decomposed
as follows This solution does not
explicitly include the control inputs or boundary condition information into
the governing equation. In order to do so an approximation of the partial
derivatives is carried out including the control inputs and the boundary data.
If *h* denotes the step size between
the points on the uniform Cartesian grid used for the finite-difference
solution, then we have [10] After substitutions can be approximated as a linear combination of
POD modes when the 's are solved in the following system model.
Then, the temporal model for the system is given by [10] where and the matrices *A* is *B* is *N* and *F* are both vectors The output equation will be simply chosen to be In this model, the dimension of the
state vector is 40 which corresponds to 40 POD modes. The first 8 POD modes
corresponding to the first 8 temporal coefficients are shown in Figure 3. The
first model corresponds to the baseline mode and the remaining modes to
actuated modes, that is, the modes due only to the input.

To assess the validity of the POD model the following test inputs, which are different from the inputs used to generate the snapshots, are applied at the boundary In Figure 4, dashed lines denote the linear combination of POD modes restricted to the boundary. Solid lines denote the boundary test inputs. As can be seen in Figure 4, there is very good agreement between the boundary conditions specified for the full-order system and the linear combination of POD modes restricted to the boundary.

The goal of model reduction is to construct another nonlinear system [10, 11] where and such that the behavior of the two systems is similar for states in some region of the state space. The reduced model is derived via the construction of an immersion/projection pair [21] where is the identity matrix, resulting in the following reduced model: This is carried out by developing an empirical balanced truncation algorithm which is based on simulation input-output measurements of the nonlinear Galerkin model. To do so we need first to introduce balanced truncation model reduction for linear time-invariant systems.

#### 3. Eigensystem Realization Algorithm (ERA)

Several frequency domain identification techniques are used in practice to identify the model parameters. One such method is the eigensystem realization algorithm (ERA) technique [26]. The ERA-based system realization model is created directly from empirical data and frequency domain characteristics of transfer functions. This method is applied to discrete time versions of system models.

A basic relationship between the Markov parameters and the input and output relationship in discrete time is An alternative form to (29) can be created not using the actual outputs and inputs but replacing the output term by the cross-correlation between the inputs and the corresponding outputs [26] where the length of the data sequence is The basic process for finding the Markov parameters starts using the ratio of the power spectral density of the cross-correlation between the inputs and outputs and the power spectral density of the autocorrelation between the input signals. These power spectral densities are given by the following: The ratio of the two power spectral densities is the frequency response function and is given by [26] Then, the final step is to take the inverse Fourier transform to find the pulse response (Markov parameter) matrices [26] The Hankel matrix containing the Markov parameters is of the following form: The individual 's correspond to the following sequence: In some cases, the input data for the ERA method might be provided by an experiment on a real system. However, in this paper a unique approach of using the Galerkin model in the place of the real system was used to generate the empirical data. The full-order system model was created using finite-difference methods. Recall that the control inputs were explicitly placed in the boundary conditions because the control inputs do not show up explicitly in two-dimensional Burgers' equation. However, the weak Galerkin model results in a nonlinear state space model that simplifies the relationship between the input and outputs. The chirp signals used for the excitation of the Galerkin model are of the following form and are shown in Figure 5:

#### 4. Empirical Balanced Truncation

The dynamics of (finite dimensional) linear time-invariant (LTI) systems are governed by a state space model of the form where is the state vector, is the input, and is the output.

The first step in applying balanced truncation is to compute a coordinate transformation such that the controllability and observability gramians, denoted and respectively, of the system are equal and diagonal. Assuming that (39) is stable, these gramians are solutions of the following Lyapunov equations [16, 17, 19]: where denotes complex conjugate transpose.

A balanced realization needs a
similarity transformation *M* such that
the transformed gramians are equal and satisfy [16, 17, 19] where the matrix is a diagonal matrix containing constants in
monotonically decreasing order.

Balanced model
reduction requires the knowledge of the controllability and observability
gramians. The latter are obtained by solving Lyapunov equations, which is
prohibitive for large-scale systems. For a system with *n* states, the controllability and observability matrices are symmetric matrices and therefore solving each
one of them involves
finding unknowns. An
alternative is to develop a balanced truncation algorithm based on empirical
gramians, which are constructed solely from a single simulation using a sufficiently
rich input. In the next paragraph we discuss how this is carried out.

The
computation cost to solve large Lyapunov equations for the controllability and observability gramians prompts us to propose a
balanced truncation algorithm, based on empirical gramians constructed from
input-output data measurements. To this end, let us first introduce the -step observability and *q*-step controllability matrices [15] which
give rise to the -step observability
and *q*-step controllability gramians As
the numbers *q* and approach infinity, these empirical gramians approach the true
gramians [15] The
goal is to find a balancing transformation matrix *M* that will approximately balance the empirical gramians, that is, The
matrix *M* can then be applied back to
the original system model to produce an approximately balanced realization.

The product of the *l*-step controllability and the *q*-step observability matrices gives a
Hankel matrix, simply a matrix that has the *i*th column identical to the *i*th row, denoted
containing the Markov parameters
of the system in the following way: for
integers *l* and *q* chosen such that [15] In
terms of the SVD decomposition of The
balancing transformation *M* is
constructed as [4] A
straightforward computation shows Balanced
truncation can be realized by
the usual way; if for some *r*, then we can partition as where A
columwise conformal partition of and yields
the immersion/projection pair [4] and
from which a reduced-order *r*-dimensional
model with state matrices is deduced The
above construction only requires estimates of the Markov parameters The
Markov parameters can be computed from a single simulation in which a
sufficiently rich input signal is applied and the output responses are
collected. In the next section, the discrete Fourier transform (DFT) is used to
map time domain data into spectral densities from which frequency response
estimates are calculated using the ERA [26].

#### 5. Application to the Galerkin Model

The empirical balanced truncation based on linear systems is applied to the Galerkin model which has an equilibrium in steady state, denoted by The rationale for doing so is that linear subspace approximations to exact submanifolds associated with nonlinear controllability and observability require only standard matrix manipulations utilizing simulation/experimental data as explained in [4, 20, 21]. The computational advantages of the scheme presented here carry over directly to the nonlinear setting.

The reduced-order model is derived as discussed through the construction of the immersion/projection nonlinear system pair This results in the following reduced-order model: If (57) has a linearization around the steady state equilibrium where and the reduced system linearization around the steady state equilibrium where The linearization of both models about equilibrium is related by Empirical balanced truncation applied to the 40th-order Galerkin model resulted in 14th-order reduced models. The first 8 temporal coefficients of the 14th-order reduced model and 2000th full-order model are plotted in Figure 6. Figure 6 shows good agreement between the temporal coefficients.

In Figure 7, we compare the Hankel singular values of the 2000th full-order linearized and reduced 14th-order empirical models. As expected the Hankel singular values corresponding to the reduced-order model are smaller than the full-order model; nevertheless, the figure shows that they are close.

In Figure 8, we compare the full-order solution of Burgers’ equation with the solution based on the 14th-order ERA model The figure shows that they behave similarly especially at the boundary where control is applied.

#### 6. Controller Design

An controller was designed based on the
linearized 14th reduced model and applied to the full-order model using Matlab.
The performance was to achieve
tracking a fixed reference signal specified for the full-order model. The
tracking problem is depicted in Figure 9, where *C* is the controller and *P* the
plant. The computation of the is based on the 14th-order reduced model.

From Figure 9, for tracking purposes the controlled output is chosen to be the
error signal *e* which is defined to be
the difference between the reference and the actual output that is, The objective of the controller *C* is to stabilize the closed-loop system and minimize the error *e*. From Figure 9, in terms of transfer function matrices of *P* and *C*, the transfer
matrix from to *e* is given by the sensitivity function defined by The control design reduces to the following
optimization. Find *C* such that the closed-loop system is stable and
minimizes Projecting onto the POD basis yields the tracking
coefficients for the reduced-order model. After computing the controller we close the loop on the original 2000th
full-order model. The controller is only 14th order since based on the 14th-order
reduced model. The projected reference onto
the POD basis initial condition and reference are shown in Figure 10. The
controlled flow with the action of the boundary controller is shown in Figure
11. The figure shows good tracking performance.

#### 7. Conclusion

Empirical balanced truncation has been considered in conjunction with POD as an approach for deriving reduced-order models and applied to 2D Burgers’ equation. Like POD, empirical balanced truncation is based on simulation/experimental data and can be implemented via standard matrix computations. Improvements to the scheme originally proposed in [4, 20] have been presented that lead to reduced data requirements that may become significant for applications such as aerodynamic flow control. Essentially, the balancing transformation is constructed via Markov parameters that can be identified from measurements collected in a single experiment/simulation. The approach has been applied with favorable results to 2D Burgers’ equation, a partial differential equation in two spatial dimensions that possesses features comparable to the Navier-Stokes equations governing fluid flow. An feedback flow controller was designed based on the empirical reduced model to achieve flow tracking. The closed loop on the full-order model shows good flow tracking performance.

#### Acknowledgments

This work was supported in part by an Air Force Summer Faculty Fellowship Program, Air Force Office of Scientific Research under Contract AF-FA9550-08-1-0450, and National Science Foundation under Contract NSF-CMMI-0825921.

#### References

- K. Willcox and J. Peraire, “Balanced model reduction via the proper orthogonal decomposition,”
*AIAA Journal*, vol. 40, no. 11, pp. 2323–2330, 2002. View at Publisher · View at Google Scholar - S. M. Djouadi, R. C. Camphouse, and J. H. Myatt, “Optimal order reduction for the the two-dimensional Burgers' equation,” in
*Proceedings of the 46th IEEE Conference on Decision and Control (CDC '07)*, pp. 3507–3512, New Orleans, La, USA, December 2007. View at Publisher · View at Google Scholar - M. Ilak and C. W. Rowley, “Reduced-order modeling of channel flow using traveling POD and balanced POD,” in
*Proceedings of the 3rd AIAA Flow Control Conference*, vol. 2, pp. 1–11, San Francisco, Calif, USA, June 2006. - D. A. Lawrence, “Empirical model reduction for active closed-loop flow control,”
*Final Report*, Faculty Summer Research Opportunity, ARL, State College, Pa, USA, 2004. View at Google Scholar - F. Leibfritx and S. Volkwein, “Numerical feedback controller design for PDE systems using model reduction: techniques and case studies,” in
*Real-Time PDE-Constrained Optimization*, SIAM, Philadelphia, Pa, USA, 2006. View at Google Scholar - C. W. Rowley, T. Colonius, and R. M. Murray, “Model reduction for compressible flows using POD and Galerkin projection,”
*Physica D*, vol. 189, no. 1-2, pp. 115–129, 2004. View at Publisher · View at Google Scholar - M. N. Glauser, M. J. Young, H. Higuchi, C. E. Tinney, and H. Carlson, “POD based experimental flow control on a NACA-4412 airfoil,” in
*Proceedings of the 42nd AIAA Aerospace Sciences Meeting and Exhibit*, pp. 4861–4870, Reno, Nev, USA, January 2004, paper no. 2004-0575. - H. A. Carlson, M. Glauser, and R. Roveda, “Models for controlling airfoil lift and drag,” in
*Proceedings of the 42nd AIAA Aerospace Sciences Meeting and Exhibit*, pp. 5627–5638, Reno, Nev, USA, January 2004, paper no. 2004-0579. - K. Cohen, S. Siegel, T. McLaughlin, and J. Myatt, “Proper orthogonal decomposition modeling of a controlled Ginzburg-Landau cylinder wake model,” in
*Proceedings of the 21st International Communications Satellite Systems Conference and Exhibit*, Reno, Nev, USA, January 2003, paper no. 2003-2405. - R. C. Camphouse and J. H. Myactt, “Reduced order modelling and boundary feedback control of nonlinear convection,” in
*Proceedings of the AIAA Guidance, Navigation, and Control Conference*, vol. 1, pp. 289–297, San Francisco, Calif, USA, August 2005, paper no. 2005-5844. - R. C. Camphouse, S. M. Djouadi, and J. H. Myatt, “Feedback control for aerodynamics,” in
*Proceedings of the 45th IEEE Conference on Decision and Control (CDC '06)*, pp. 1–14, San-Diego, Calif, USA, December 2006. - W. R. Graham, J. Peraire, and K. Y. Tang, “Optimal control of vortex shedding using low-order models—part I: open-loop model development,”
*International Journal for Numerical Methods in Engineering*, vol. 44, no. 7, pp. 945–972, 1999. View at Publisher · View at Google Scholar - H. Banks, R. Del Rosario, and R. Smith, “Reduced order model feedback control design: numerical implantation in a thin shell model,” Tech. Rep. CRSC-TR98-27, Center for Research in Scientific Computation, North Carolina State University, Raleigh, NC, USA, June 1998. View at Google Scholar
- E. Caraballo, M. Samimy, and J. DeBonis, “Low dimensional modeling of flow for closed-loop flow control,” in
*Proceedings of the AIAA Aerospace Sciences Meeting*, Reno, Nev, USA, January 2003, paper no. 2003-0059. - W. J. Rugh,
*Linear System Theory*, Prentice-Hall, Upper Saddle River, NJ, USA, 2nd edition, 1996. - K. Zhou, J. C. Doyle, and K. Glover,
*Robust and Optimal Control*, Prentice-Hall, Upper Saddle River, NJ, USA, 1996. - K. Zhou and J. C. Doyle,
*Essentials of Robust Control*, Prentice-Hall, Upper Saddle River, NJ, USA, 1998. - A. Defant and K. Floret,
*Tensor Norms and Operator Ideals*, Elsevier Science, Amsterdam, The Netherlands, 1993. - B. C. Moore, “Principal component analysis in linear systems: controllability, observability, and model reduction,”
*IEEE Transactions on Automatic Control*, vol. 26, no. 1, pp. 17–32, 1981. View at Publisher · View at Google Scholar - J. Hahn and T. F. Edgar, “Reduction of nonlinear models using balancing of empirical gramians and Galerkin projections,” in
*Proceedings of the American Control Conference*, vol. 4, pp. 2864–2868, Chicago, Ill, USA, June 2000. - M. Condon and R. Ivanov, “Empirical balanced truncation of nonlinear systems,”
*Journal of Nonlinear Science*, vol. 14, no. 5, pp. 405–414, 2004. View at Publisher · View at Google Scholar - C. A. Antoulas,
*Approximation of Large-Scale Dynamical Systems*, Advances in Design and Control, SIAM, Philadelphia, Pa, USA, 2006. - V. Sreeram, “Frequency response error bounds for time-weighted balanced truncation,” in
*Proceedings of the 41st IEEE Conference on Decision and Control*, vol. 3, pp. 3330–3331, Las Vegas, Nev, USA, December 2002. - T. J. Chung,
*Computational Fluid Dynamics*, Cambridge University Press, Cambridge, UK, 2002. - P. Holmes, J. L. Lumley, and G. Berkooz,
*Turbulence, Coherent Structures, Dynamical Systems and Symmetry*, Cambridge University Press, Cambridge, UK, 1996. - J. Juang,
*Applied System Identification*, Prentice-Hall, Upper Saddle River, NJ, USA, 1994.