An optimal input design procedure

Abstract

This technical communique presents a method which aims at improving input design searching in some dynamical systems. The original idea of the proposed method is the combination of a dynamic programming method working on square wave inputs followed by a Quasi-Newton algorithm working on infinite differential switching mode inputs. These infinite differential functions approximate square wave inputs in order to enlarge the set of admissible inputs while ensuring a reasonable computation requirement which is not the case of classical methods based on dynamic programming only. Moreover, they correspond to practical inputs. The precise description of the approach is followed by an application in aerospace sciences.

Introduction

A large part of the literature is devoted to optimal input design in industrial systems (in marine systems, for example, Fossen, 2002 or in aircraft systems, for example, Chen, 1975, Mehra, 1974, Morelli, 1999). The nonlinear controlled dynamic systems considered in this paper can be rewritten as

$$\dot{x}(t,p)=f(x(t,p),p)+u(t)g(x(t,p),p)\text{,}$$

$$y(t,p)=h(x(t,p),p),x(0)=x_0\text{,}$$

where $x(t,p)\in {\mathrm{R}}^n$ and $y(t,p)\in {\mathrm{R}}^m$ denote, respectively, the state variables and the measured outputs. The input function u is assumed to be piecewise continuous or differentiable. The vector $p=(p_1,\dots ,p_l)$ represents the parameters to be estimated and belongs to an admissible set $\mathrm{\Omega }$, subset of ${\mathrm{R}}^l$. The time is assumed to belong to $[0,t_{\max }]$. The functions $f(x,p),g(x,p)$ and $h(x,p)$ are real and analytic on M for every $p\in \mathrm{\Omega }$ (M is a connected open subset of ${\mathrm{R}}^n$ such that $x(t,p)\in M$ for every $p\in \mathrm{\Omega }$ and every $t\in [0,t_{\max }]$). The single-input case is considered for notational simplicity; all the results can be readily generalized.

The input u is designed by considering an approach based on two successive steps.

The first uses dynamic programming for minimizing a cost function linked to the Fisher information matrix. In this step the admissible set of inputs is built with full amplitude square waves. It leads to a first input $\widehat{u}$. A global exhaustive search can be carried out but the computational requirement, specially computation time, expands rapidly by increasing the choice of square wave input and requires a restrictive admissible square wave input set.

In order to enlarge the admissible set, the square wave inputs are approximated in $L_2([0,t_{\max }])$ by differential switching mode inputs, which involves some new parameters. Thus the cost function becomes regular and the second step consists in minimizing it by a Quasi-Newton algorithm starting from an approximation of the result $\widehat{u}$ given by the first step. By this way all the admissible inputs can be taken into account. A description of this method is presented in the following section.

The measurement data are assumed to be given by

$$y_m(t_i)=y(t_i,p)+\nu (t_i),i=1,\dots ,N\text{.}$$

The measurement noise $\nu (t_i)$ is assumed white gaussian with zero mean and $E[\nu (t_i){\nu }^{\top }(t_j)]=R{\delta }_{\mathit{ij}}$ where R is the measurement noise covariance matrix and $i,j=1,\dots ,N$.

The test time T is assumed to be fixed and such that $T⩽t_{\max }$. The integer N is the total number of sample times.

Our purpose will be illustrated by an example concerning the longitudinal motion of a glider. The projection of the general equations of motion onto the aerodynamic reference frame of the aircraft and the linearization of aerodynamic coefficients (Wanner, 1984) give the following system:

$$\left\{\begin{array}{c}\dot{V}=-g\sin (\theta -\alpha )-{\displaystyle \frac{1}{2m}}\rho {\mathit{SV}}^2(C_x^0+C_{x\alpha }(\alpha -{\alpha }_0)\hfill \\ +C_{x{\delta }_m}({\delta }_m-{\delta }_{m_0})),\hfill \\ \dot{\alpha }={\displaystyle \frac{2}{2\mathit{mV}+\rho {\mathit{SlVC}}_{z\dot{\alpha }}}}\left\{\mathit{mVq}+\mathit{mg}{\displaystyle \frac{\cos (\theta -\alpha )}{V}}\right.\hfill \\ -{\displaystyle \frac{1}{2}}\rho {\mathit{SV}}^2\left(C_z^0+C_{z\alpha }(\alpha -{\alpha }_0)\right.\hfill \\ \left.\left.+C_{\mathit{zq}}{\displaystyle \frac{\mathit{ql}}{V}}+C_{z{\delta }_m}({\delta }_m-{\delta }_{m_0})\right)\right\},\hfill \\ \dot{q}={\displaystyle \frac{1}{2B}}\rho {\mathit{SlV}}^2\left\{C_m^0+C_{m\alpha }(\alpha -{\alpha }_0)+C_{\mathit{mq}}{\displaystyle \frac{\mathit{ql}}{V}}\right.\hfill \\ +C_{m\dot{\alpha }}{\displaystyle \frac{2l}{2{\mathit{mV}}^2+\rho {\mathit{SlV}}^2{C}_{z\dot{\alpha }}}}\left[\mathit{mVq}\right.\hfill \\ +\mathit{mg}{\displaystyle \frac{\cos (\theta -\alpha )}{V}}-{\displaystyle \frac{1}{2}}\rho {\mathit{SV}}^2\left(C_z^0+C_{z\alpha }(\alpha -{\alpha }_0)\right.\hfill \\ \left.\left.+C_{\mathit{zq}}{\displaystyle \frac{\mathit{ql}}{V}}+C_{z{\delta }_m}({\delta }_m-{\delta }_{m_0})\right)\right]\hfill \\ \left.+C_{m{\delta }_m}({\delta }_m-{\delta }_{m_0})\right\},\hfill \\ \dot{\theta }=q.\hfill \end{array}\right.$$

In these equations, the state vector x is given by $(V,\alpha ,q,\theta {)}^{\top }$, the observation y is full ($y=x$), the input u is ${\delta }_m$ ($u_0={\delta }_{m_0}$) and p the parameter vector to be identified is equal to the dynamic stability derivatives $(C_{z\dot{\alpha }},C_{\mathit{zq}},C_{m\dot{\alpha }},C_{\mathit{mq}})$. The real V denotes the speed of the aircraft, $\alpha$ the angle of attack, ${\alpha }_0$ the trim value of $\alpha$, $\theta$ the pitch angle, q the pitch rate, ${\delta }_m$ the elevator deflection angle, $\rho$ the air density, g the acceleration due to gravity, l a reference length and S the area of a reference surface. Experiments are based on free flights of scaled models in a laboratory.

In the case of aerospace models, constraints arising from practical flight test considerations were imposed on all input amplitudes and selected output amplitudes. Control input amplitudes are limited by mechanical stops, flight control software limiters, or linear control effectiveness. Then it is given by

$$|u(t)-u_0|⩽\mu \forall t\text{,}$$

where $\mu$ is a positive constant, and $u_0$ the trim value of u.

Selected output amplitudes must be limited to avoid departure from the desired flight test condition and to ensure validity of the model. In addition, constraints may be required on aircraft attitude angles for flight test operational considerations, such as flight safety.

This paper is organized as follows. In Section 2, the approach based on two successive steps is given in order to calculate optimal input functions. Then it is applied to the input design of the aircraft model (3), which is presented in Section 3.