- Email: [email protected]
State Synthesiser; A Digital Observer for Spacecraft Attitude Control Systems
Copvright
© IFAC Automatic Contro l in Space
~oo rdwijke rhout.
The
~rtherlands
1982
STATE SYNTHESISER; A DIGITAL OBSERVER FOR SPACECRAFT ATTITUDE CONTROL SYSTEMS S.
J. Dodds
Attitude Control Systems Group, Satellite Ditq·sioll , Marconi Space and Defence Systems Ltd., Browns Lane, Th e Airport , Portsmouth , Hampshire , UK
Abstract. The state synthesiser is a form of digital state observer which is optimised for speed of response rather than minimum noise content. The serious limitation often imposed by a Kalman filter on the speed of response of a high preCision control loop to transient disturbance torques (such as produced by orbit change manoeuvres) is removed by switching from the Kalman filter output to the state synthesiser output when the pointing error exceeds a given threshold. In the state synthesiser, the state transition equation for the dynamics is repeatedly applied to construct an estimate of the present state given a finite, minimum length, sequence of past measurements together with the control function operating during the sequence. Simulations are given of a control loop embodying a statE' synthesiser, steady state Kalman filter and sub-time optimal control law (which is adaptive to a disturbing torque estimate) for attitude control of a spacecraft with a flexible structure. A lumped parameter dynamics model is utilized. The principle is not limited to the above application and may be relevant to other time invariant, linear plants of arbitrarily high order. Keywords. Adaptive control; Aerospace computer control; A tti tude control; Bangbang control; Digital computer applications; Discrete systems; Observers; State estimation. perfect performance is obtained in absence of plant and measurement noise or parameter mismatching.
INTRODUCTION The principle function of the state synthesiser is to provide a state estimate which follows, as closely as possible, a rapid Change in the plant (dynamics) state produced by an unexpected disturbance of unknown magnitude. This enables an automatic control loop to respond rapidly, reducing transient errors.
The initial motivation for development of the state synthesiser is the need to minimise transient pointing errors in high precision spacec raft, due to disturbing torques occ urring during firing of orbit adjustment thrusters. Such a problem is encountered on the X-ray astronomy satellite, Exosat. The latter spacecraft employs an alternative scheme embodying gain switching of a more conventional state estimator as described by Dodds (1981 a). Recognition of the dangers of impairing the stability of high order control loops by employing such gain switching, however, leads to the concept of running two state estimators concurrently, one with stochastically optimised gains for fine pointing and the other with higher gains set for fast response. A sudden disturbance causes the optimal state estimator error (the difference between the measurements and corresponding quantity derived from the estimator dynamics model) to grow rapidly at first. The passing of the error beyond a given magnitude threshold is used to indicate the onset of the disturbance and the state input to the control law is switched from the optimal estimator to
The concept is eminently suitable for application where the control law requires an estimate of the complete state vector. Additional state variables which characterize the disturbance are also estimated, as a necessary part of the process of obtaining the plant state estimate and these may also be utilised by the control laws. The controlled plant is assumed time invariant, or, at most slowly time varying but may be of arbitrarily high order, have multiple control inputs and no restrictions on open loop pole locations. The state synthesiser functions in a similar fashion to a finite impulse response filter and employs a non-recursive algorithm in which the estimate of each state variable is obtained as a weighted sum of a finite sequence of past measurements and corresponding control signals. As with the Kalman filter, mathematically
A: ':
62 _ :' "
191
s. J.
192
the 'fas t' e s ti ma to r as lo ng a s the stat e e s tim ato r e rr o r magni t ude rema in s a bove t he thr esho ld, enab li ng con tr o l t o rqu es t o be prod uced rapid ly to combat th e d is t urba nce . Such a sc heme is il lus tr ated in F ig . I. Inc rease o f t he speed of respo nse o f a d ig ita ll y imp leme nted, disc re t e es ti ma to r by inc rease o f ga in s is u lti mately li m ited by t he it e ra t io n pe ri od , h, uns t able osc ill a t io ns o f th e s t a t e est ima t e occu rr ing a t a f r e qu e ncy o f 1/2h Hz if t he ga in lies outside a ce rta in st ab ilit y boundary : This is due to closed loo p po les o f the es t ima to r becoming compa rable in magni t ude wi t h I /h, as th e gains a re inc reased . The sta t e svnthesise r overcomes th is limi t a t ion and also enab les an exact state es t imate to be obtained with a r bitra r y initia l conditions in t he minimum number of algo r ithm it e ratio ns.
Dodds The disturban c e, V( t), is assumed to be deterministi c in na ture and c omposed of c o m ponents whi c h can be ap pro ximate d by so lutions to the st a t e equat io n
W = QW
(3 )
whe re W is a st a t e vec t o r, o f d im e nsio n, w, a ssoc iated with th e d is turb a nce and Q is a fi xe d matrix c har ac teri z ing the di s turban ce process . E xa m ples inc lude st e p di s turbing t o rques and cyc li c d istu rbing t o rqu e s appr ox imate d by a t ru nca t e d Fou r ie r se ri es . The d ist urba nce vec t o r is then ( 4)
V = PW
where P is a fix e d ma tri x o f di me nsio n v x w.
,A, fas t s t a t e observe r , such as th e state svnthes ise r, is necessarily se nsi ti ve to measu remen t and plant noise a nd is t he re fo re int ended to supp lement an op ti mal st a t e estima to r (Kalman fil te rl. In si t ua t io ns where the measu reme nt and plan t noises a re no t too seve re , however , the s ta te synt hesise r may ent ire ly replace the Ka lman fil te r.
St oc has ti c d is t urbances are no t inc lud e d in eq ua ti o ns (1) a nd (2), since t he de rivati on o f th e st a t e s ynthes ise r a lgo r ithm is de t e rmin is ti c in na t ur e and does no t re qu ire kn o wledge o f plan t and me a sur e me nt no ise. The latte r a re inc luded in the simu la ti ons , however, in o rd e r t o assess se ns iti vi t y t o rea li s t ic senso r a nd co ntr o l ac tu a t o r noise .
The simu la t ion results presented in t his pape r a re fo r a tt i tude control about a s ingle axis o f a spacec raft embodying a sing le , domina nt , f lexu re mode. This re lat ive ly s imp le examp le serves to demonst r ate t he ope r at io n o f th e s t a te sy nt hes ise r , but it mus t be rea li sed th a t t he proced ur es given fo r ca lc ul a t io n o f t he s t ate sy nt hesiser matrices a re ge ne r a ll y applicab le and may inc lude simul t aneous t h reeaxis control with linear inte r- ax is co up ling, many f lexu re modes a nd m ulti ple con tr ol ac tua to rs pe r ax is .
Fi na lly , e qu a ti ons (1), (2 ), (3 ) and (4 ) may be comb ined so that , in pa rt it io ned fo r m ,
U
+
(M OJ
Y
+
DU (5 a )
a nd in the sta nd ard fo r m, DEFI NITI O"-! OF PLAN T AND DI STUR BANC E PARAME T ERS In view o f the app licab ili ty t o a ny ti me in va ri a nt, linea r plan t, th e concep t is deve lo ped in ge ne ra l t e r ms , app li ca ti on t o a spec ifi c spacec ra ft dy nam ics being au t oma ti ca ll y included in this general fo r mula t ion. The plant to which the state sy nthesise r app li es !s assumed to be of t he fo r m
Z = FZ
{
x
AX
+ BU
Y
ex
+
(5b) DU
where X is an ove ra ll sta t e vec t o r of d im e ns ion N = n + wa nd A, B a nd C ar e pla nt, inp ut a nd measu reme nt ma tr ic es co rres ponding t o th e pa rtiti oned ma tri ces of e quation (5a). DEVElOPME NT OF THE C O NC EPT
+
GU
+
EV
whe re Z is th e pla nt st ate vec t o r , of di me ns io n n, U is th e co ntr o l vecto r o f d imen sio n r, V is th e di s turban c e vec t o r of d im e ns io n v, F is the plant matri x, B is th e input matrix and E is t he distur banc e ma tr ix. Th e me a sur e ment equ a t io n is of the fo r m Y = MZ + DU
( 2)
where Y is the measure ment vector of dimension m, M is the outp u t matrix and D is the feed-forward matrix.
In o rd e r for the c ur re nt s tate to be estimated in minimum time, witho ut pr ior kno wled ge of the state, the minimum quantit y of input / output obser vation s of the plant must be ut ilized. In the state synthesiser, these observations take the form of a stored past histor y of y( t) and U( t) , sampled at the algorithm iterat ion interval, h. The essential fa c t which determines the minimum number of samples of (Y, U) is that a sufficiently large set of temporariy spaced measurement coordinates serve as well as spacially separated co-ordinates within . Y for the purpose of re-
State Synthesiser constructing the state vector at a given time. In fact, as proven in the algorithm derivation below, the product of the number of samples of Y and the dimension of Y must be at least equal to the order of the system described by equation (5).
193
In general, C is of dimension m x N and is non-invertible, but it is possible to construct an invertible measurement equation by utilizing more than one sample of Y and relating them to X k by repeated application of the state transition equation, as shown below. Suppose that s consecutive input/output samples are utilized, including the current sample (Y k' lI k _I )·
DERIVA TlON OF THE STATE SYNTHESISER ALGORITHM
Then This present derivation assumes constant values of U between algorithm updates. Generalisation to known U(t) between algorithm upda tes, however, should be possible, but is not treat e d in this paper. The state transition equation and measurement equation corresponding to the continuous state equations (5) is Xk '"
(6a)
Y k '" CX k + DU k _ 1
(6b)
Y k-q
=
f
h
o
+
(9)
DlIk_q_1
and inverting the state transition equation: X _ '" ~-I k q
Xk_q+ I _~-I
'V
lIk_q
'" ~-2
X k _q +2 - ~ -2
'f
lIk_q+ I
- ~ -I
'f
(10)
lIk_q
Repeated application of this substitution yields
where cP is the state transition matrix and 'V will be termed the input transition matrix, so that 'I' (h)
CX k _q
)('_q '" "K
J,
'!'
-q X
k
q L: ,;, - i .\'t', I Ik-q + i-I 'I' i =I
-
(11 )
q '" 0, I, 2, .•... S-I
(h-'t ) Bd't.
Substituting for Xk q equation (11) Yields -
The state synthesiser algorithm is derived simpl y by repeated application of equation (6) to past samples of (Y, U).
in
equation
(9),
using
q C ~-q ~
Now, it is assumed, as is usually the case, that m < N where m is the dimension of Y and N is the dimension of X. for m '" N, the state synthesiser algorithm is trivial and consists of the inverted measurement equation
Z i=I
+
D lIk _q _ 1 (12)
q '" 0, I, 2, ..... S-I
(8)
(sm x I) ~
Yk Yk - I
Now, equation (12) may be written for every sample, (y k-q' lI k _q _ I ), in partitioned matrix form as follows:
(sm x N) C
----
X k
C~-I
-D
CqJ-I \11
0
0
-D
0
~
.............. , "-
Yk - 2
C~-2
I I
I I I
I t
I
I
I
Yk - s + 1
C~-2
I
I
C~
-(s-I)
(sr x I)
(sm x s r)
~
I
If'
Cq,-I'ji
,
"-
"-
-D
'-"-
"-
cf(s-l)~ Cf(S-2)~
"-
"-
"-
"-
"-
"-
"-
..•..•. '
"-
"-
'--
"-
"-
"-
0
lI _ k 1
0
lIk_2
,
,
0 "-
Cfl~
"-
"-
-D
Uk _ s
(I3)
s. J.
194
The expressions written in parenthesis over equation (13) indicate the matrix dimensions. Now, equations (13) represent s.m linear simultaneous equations in the N unknown components of Xk , and can be solved as they stand only if s.m = N. If s.m < N, there are insufficient equations, tantamount to the (Y) samples containing insufficient information for reconstruction of the state vector. Hence the minimum number of samples can be expressed as
Dodds amount equal to the time taken by the control law to calculate and issue U to the control k actuators. This problem may be avoided by ,A
applying the state transition equation to X in k order to predict the state one iteration step ahead in time. Thus
(18) (14)
smin = int(N/m)
where int(x) is defined as the smallest integer which is ~ x. [n many cases, N/m is not an integer, in which case equations (13) are overdetermined for s = smin. This problem may be overcome, however, by lntroduclng a weighting matrix, W, of dimension N x s.m. First, equations (13) are written in compact form as
I
I
I
issue
of
the
latter
to
the
control
actuators being delayed until sampling instant k+ I. may
be
re-written
in compact
Where [UJ k has been re-defined to include U as its first
r elements.
Sand S
y
u
k are the
fixed, predicting state synthesiser matrices of dimension N x s.m and N x (s+ 1).r, respectively. In terms of the matrices of equation (18),
(Nxsr)
T'I) TI) I
(19)
(15)
where the matrix dimensions are indicated in parentheses. Both sides of equation (15) are then pre-multiplied by W to yield (NxN)
Uk+ I,
~ srxI )
[yJ k = [C]Xk + [DJ[UJk
(Nxsm)
step is allowed for computation of Xk+ I and
Equation (18) form as
( smx s r ) ( smx I ) ( smxN)
Equation (18) constitutes the algorithm of the predicting state synthesiser which is the version recommended for use in automatic control applications. In this case, one iteration
( 1" ' 1 )
(20)
w[yJ k= [W[CJjXk + [WeD] ][U]k
(16)
in which the matrix, W [C] , may be inverted. [t can br shown that for observability of X, Rank [C J = Nand it follows that W may be chosen so that W may be inverted.
Cc J
The basic state synthesiser algorithm is then obtained from equation (16). The synthesised A
state is denoted X in order to distinguish it /\
from X, the symbol usually used to denote the state estimate from a Kalman filter. Then, prfmultiplying both sides of equation (16) by
(21) In absence of plant parameter mismatch and stochastic disturbances, the predicting state synthesiser is capable of observing the state perfectly, after s+ I initial iterations, occupied by an acquisition transient, analogous to the behaviour of a finite impulse response filter. The algorithm is particularly simple in that it is non-recursive and the estimate of each state variable is obtained independently as a weighted sum of a minimum length sequence of past measurements and control values.
LW[CJ] -I yields FIL TERING STATE SYNTHESISER
For
the
cases
in
which
smin
<
N/m
(see
equation (14)), there are m.smin - N redundant ~
Now, Xk is the estimate of the state at the instant of the current measurement sample, y k' and since the next control, Uk, is calculated by the control law as a function of Xk , the control would be applied late by an
measurement coordinates. Suitable choice of W then enables the redundant measurements to be used to advantage, giving the state synthesiser filtering properties. In fact, there is no theoretical upper bound on the number of samples, s, of (Y, U) and in cases where stochastic performance of the minimum
State Synthesiser sequence
(s
smin)
state
synthesiser
is
unsatisfactory, a trade off may be made by increasing s at the expense of increased acquisi tion time «smin+ I)h increased to (s+ I)h). It may further be shown that W may be chosen so that for any s < N / m, the state synthesiser becomes a least squares estimator in the following sense: Let Yk(t) = [C]xk(t) where tk (s+ I)h ~ t ~ tk and Xk(tk) ~ Xk, then if (22)
then Y k(t) is a least squares fit through the measured values Y i' i = k-s, k-s+ I, ..... k-l, k, as
generated
assuming the k values U i , i=k-s+l, k-s+2, ..... k-I, k.
xl x2
from
I:
X
0
0
x3
0
0
0
x4
0
0
_ w
Xs
0
0
0
and
z1
w
control
0
0
0
0
0
2
0 0
0
2
SPACECRAFT ATTITUDE CONTROL APPLICATION For the purpose of demonstrating the state synthesiser, a simple spacecraft example is considered, embodying a single flexure mode. The simulation is confined to gas-jet control, about a single axis, of a spacecraft comprising a main rigid body and a significantly flexible and rigidly attached appendage, such as a solar array, the whole craft exhibiting a single dominant vibrational mode. The flexible mode is represented by a second rigid body with one degree of rotational freedom relative t o the main body, the two bodies being restrained about the control axis by means of a linear torsion spring with zero damping. The gas jets are assumed to apply control torques direc tly to the main rigid body. The relevant plant equation, in the companion form, is
"1 x2
0
x31
0
J
0
w
2
0
l+ A 0
z3
0
-I
0
1I
lo
:l
xl
0 1
x3
x2
i
z4
0
0
0
-I
z sJ
0
0
0
0
y
z 1 + nm
(u + npl
1
,,J
l+A z2
0
+
0
0
19 5
:~
I
I x4
1
l xsJ
(23)
s. J.
196
where z I and z2 are
the
Dodds
pointing error and
pointing error rate, respectively, and z3 and z4 are, respectively, the equivalent flexural-mode displacement and rate. The flexural mode natural frequency is w , the control angular acceleration is u and the angle measurement is y. The moment of inertia ratio between the two bodies of the lumped parameter dynamics model is denoted by A (mode inertia to main spacecraft body inertia). The plant and measurement noises are denoted np and nm respectively.
In cases where limit cycle optimisation requires impulses which are shorter in duration than h, the state synthesiser must be developed further to accommodate U k as a known function of time, defined in the interval (t k , tk+ I). No difficulties are anticipated in this development. Where reaction wheels are used for control, however, the digital control law may calculate control torque demands intended to be held constant for each iteration, rendering the present state synthesiser directly applicable.
The deterministic disturbance is assumed to be a constant disturbing acceleration and is modelled in equations (23) as x5 and z5. Note that the state in X-space is to be estimated, since the control law utilized for the simulation is a bang-bang control law which requires, as input, the state of the plant in the companion form. This control law, which applies to any time invariant, single input, linear plant is presented in detail by Dodds (J98Ib). The state transition equation of the plant is
h
XI
( 1 - c) /w 2
si",
Digital simulation results of a spacecraft attitude control system of similar form to that depicted in Fig. I, using the dynamics equations (23), are presented in Figs. 2, 3 and 4. The computer utilizes a word length of 39 bits. Angular units are given in seconds of arc since applications to spacecraft of high pointing precision are envisaged.
h-slw)/w 2 Qh 2 +(c_1 )/w 2 )/",2
(
( 1 _ c) Iw 2
x2
0
x3
0
0
c
x4
0
0
-ws
c
0
0
0
0
x5
DIGIT AL SIMULA TlON RESULTS
slw
xI
(tIl 2 • (c - I) 1",2) 1",2
( h - si", ) I", 2
x2
(h_s/",)/w 2
(l_c)/w 2
x3
(l_c)/w 2
slw
s
x4
uk
/w
x5
k,1
k
where h is the iteration step length of the onboard digital processing for state estimation and control. .A lso, c = cos w hand s = sin who In this case, N = 5 and m = 1, so that, by equation ([4), smin 5. Now, in general, the
The nominal values of the constant parameters of the simulation are chosen as below to represent one control axis of a three axis stabilised satellite having a flexible appendage with a significant vibrational mode. Inter-axis coupling is not included.
predicting state synthesiser requires smin + 1
Control acceleration from gas jets
samples due to the last term of equation ([2), but, for this example, D = 0 so that only 5 samples of y and n are required for the
= 100 arcsec s-2
~
minimum length sequence to generate X. The predicting state synthesiser matrices, Sy and Sw are therefore, each of dimension 5 x 5 and are generated using state transition matrices similar to that used in equation (24), but with h replaced by -ih, where i = 1, 2, 3 and 4. Note also that, for this example, the control between iterations is assumed constant. A value of 120 ms is chosen for this example, this being the minimum 'on' time of the control jets.
Flexure mode natural frequency, w=12rads- l Lumped parameter model inertia ratio, 3
,,=
State synthesiser iteration interval, h = 0.125 s. The above control produced by a 0.05 operating at a 1 spacecraft moment of m 2•
acceleration would be Nm cold gas thruster m arm with a total inertia of about 400 kg
197
State Synthesiser Control loop responses Fig. 2 shows the control system responses to a disturbing acceleration step of 80 arcsec s-2 (80% of the available control torque) applied at t = 1.2 s and removed at t = 20 s.
simulation to obtain the gyro angle noise. Performance is not drastically reduced by the presence of plant and measurement noise, the peak pointing errors being of the same order as those for the ideal case of Fig. 2b. At t = 1.8 s, the disturbing acceleration estimate, C'd' jumps as the error magnitude threshold switch (rei. Fig.l) changes from the Kalman fi!1er to the state synthesiser. The variations in Lld due
For comparison purposes, Fig. 2a shows the response obtained with a steady state Kalman filter only, in which case the dynamics of the filter dominates the control loop transient, resulting in relatively large pointing errors, the peak error following application of the disturbance acceleration step being about 500 arcsec and that following removal of the disturbance being about -200 arcsec. Note, that in order to show this pointing transient on the same angle error scale as the remaining responses, in which the state synthesiser is brought into play, an unfolded cylindrical plot is utilized. The final portion of the settling transient following the -200 arcsec peak is of limited interest and has been omitted in order to accommodate the plot on the same time scale as the remaining plots.
towards zero, as the state sy nthesiser is again brought into play, following the removal of the high level disturbing acce lera tion. This is followed by a period of noisy disturbing acceleration estimate until about t = 24s, when the Kalman filter is switched in again.
Figs. 2b through 2i show control loop responses in which the error magnitude switch of Fig. I switches to the state synthesiser when the
Figs. 2d to 2i are noise free runs but with various parameter mismatc hes as desc ribed below.
Kalman filter error, magnitude of 2 arcsec.
zI
zI
exceeds
a
In Fig. 2b, the dynamics model parameters of the Kalman filter and state synthesiser are precisely matched to those of the real dynamics and no plant or measurement noise is included. The peak pointing error magnitude during application of the disturbance acceleration is reduced to about 16 arcsec, in contrast to 500 arcsec without the sta~e synthesiser. The peak pointing error following removal of the disturbance acceleration is reduced to about 5 arcsec in magnitude, in contrast to 200 arcsec without the state synthesiser. This run also serves as a standard of comparison for the simulations described below in which realistic imperfections are introduced. Fig. 2c shows the control system response with perfectly matched parameters as in Fig. 2b, but with plant noise and measurement noise. The plant noise is derived from typical cold gas thruster data and is taken as having a constant spectral density of 0.02 (arcsec s-2)2/Hz, being simulated as Gaussian
to the stochastic disturbance are evident in Fig. 2C. At approximately t = 5.5s, the state estimate from the Kalman filter is again used, as shown by the marked reduction in noise on the disturbing acceleration es timate. The end of the Kalman filter transient can just be seen, as the disturbing acce leration estimate monotonically approaches the correc t value, u . d ~ At about 21s, u d steps abruptly down
In Fig. 2d, th e con trol acceleration assumed in the state synthesiser and Kalm a n filter is overestimated by 10%. Again, at about t = 1.8s, the state synthesiser o utput is used and the initial variations of ""ud a bout t he true value, ud, are due to t he finite acquisition transient of the state synthesiser which occu rs after application of the negative contr ol jet. The difference between the control accelera t io ns applied to the re al dynamics and dynam ics model is auto matica ll y treated as a s tep disturbance acceleration of -1 00 - (-110) 10 arcsec s-2 and is added t o the estimate of ud. A
The variations in u d a fter removal of the real disturbance acceleration are due to the alternate pulsing of the over-estim ated control jets. This variation is drastically reduced when the Kalman filter estimate is used after about t = 25s. In Fig. 2e, the control acceleration assumed in the dynamics models is under-estimated by 10%. The additional disturbance acceleration due to application of the negative control jet
noise with an r.m.s. value of 0.28 arcsec s-2 over a bandwidth of 1/2h ~ 4 Hz where h is the iteration interval of 0.12 s. The gyro noise is assumed to be white in rate and, using typical high precIsion rate integrating gyro data, a somewhat pessimistic spectral density
steady state estimate, ud' is, therefore, lower
of 0.017 simulated
being r.m.s.
remarks apply to this case as those given for Fig. 2d.
value of 0.27 arcsec s-I over a bandwidth of 4 Hz. This rate noise is integrated in the
Figs. 2£ and 2g show runs in which the flexure mode natural frequency, w , assumed in the
(arcsec s-I )2/Hz is taken, as Gaussian noise with an
is then -lOO - (-90) = -10 arcsec s-2.
The
~
than ud by 10 arcsec s-2.
Otherwise, similar
s.
198
state synthesiser and Kalman filter is, respectively, over-estimated by 20% and underestimated by 5%. The control loop response is found to be insensitive to over-estimation of w but very sensitive to under-estimation of w . Runs for over-estimation of more than 20% are not presented since errors of this magnitude are unlikely to occur in practice. Increasing the under estimation by more than 5% produces relatively large peak pointing errors (about 70 arcsec at 5%) and can lead to control loop instability.
The oscillations in
~d
occur while the estimate
of the state synthesiser is being used and comparison with the rate trace, z2' shows that they occur at the flexure mode natural frequency. Comparison between F~s. 2f and 2g shows the oscillations of u to be d approximately anti-phase with the oscillations of z2 for the +2096 mismatch but in-phase for the -5% mismatch. The tendancy towards control loop instability with natural frequency under-estimation is probably associated with these phase relationships. Figs. 2h and 2i show runs in which the modelled inertia ratio, A, is, respectively, overestimated and under-estimated by 10%. These results are similar to those for the natural frequency mismatching, in that similar phase relationships occur between 'C'd and z2.
This is
because the natural frequency of the dynamics model is proportional to 1\. Note that although the -10 % m isma tch of Fig. 2i produces an acceptable peak pointing error of about 35 arcsec, there is evidence of control loop instability in the z2 trace where the Kalman filter is being used just before the removal of the disturbing acceleration step at t = 20s. State Synthesiser Acguisition Transient Fig. 3 shows the individual state variable outputs of the perfectly matched state synthesiser during application at t = 1.56s of a high level disturbance acceleration step, u . d The continuous curves are the real state variables and the 'staircase' functions represent the corresponding state synthesiser outputs. During each iteration of 0.12s duration, the predicting state synthesiser computes the encirled points on the 'staircase' traces of Fig. 3. Before application of the step disturbance acceleration, the state synthesiser produces the correct state estimates. Errors in the state estimates are evident for five iterations following application of the disturbance step, in accordance with the theory: The number of steps occupied by the acquisition transient is equal to the product of the dimensions of the measurement vector (I) and the number of state variables (5).
J. Dodds Even during the transient, the estimates tend to follow the true state variables, except for f!I.
ud
in
which
exaggerated.
the
errors
are
the
most
The unusual transient behaviour
A-
of ud' unlike that of any recursive estimator, such as the Kalman filter, is attributed to the oscillation of the elements of the associated row of the Sy matrix. This oscillation of terms is sequence
-I
,~
-2
basically
,
-3
due
•••••
to oscillations in the
which occur due to the
complex eigenvalues of Ijl. COMPARISON WITH RECURSIVE ALGORITHM With an iteration period of 0.12S, it impossible, in the flexible spacecraft to stabilise the state estimator by choice of gains, due to the imaginary
is found example, suitable poles of
the dynamics model at + 12 rad s-I. All attempts result in an oscillation, at 0.24s period of the state estimates about the true values. The problem is removed in the simulations of fig. 2 by iterating the state estimator twice for every iteration of the state synthesiser, with the dynamics model state transition matrix calculated using an iteration period of 0.06s. Comparable results to those of fig. 2b can be obtained, however, using a state estimator with increased but stochastically sub-optimal gains. This is really an unfair comparison, however, since the state estimator has the advantage of operating at twice the iteration frequency of the state synthesiser: Comparison should be made at the same interation frequency, so that in the spacecraft example simulated, the state estimator, even used as an observer with suboptimal gains adjusted for increased speed of response, proves inferior to the state synthesis er; unsatisfactory, in fact, due to the oscillations of the state estimate about the real state, as described above. A recursive algorithm, such as the steady state Kalman filter with fixed gains, involves fewer arithmetic operations per iteration than the state synthesiser but, besides the iteration period limitation referred to above, also requires a theoretically infinite number of iterations to obtain an exact state estimate commencing with arbitrary initial conditions. The state synthesiser, therefore, compare$ favourably for 'fast' state estimates such as required in the application described in this paper.
State Synthesiser Operation of the state synthesiser with extended iteration period
19 9
simulated as for the run of Fig. 2C, but with zero values of demanded control signal and zero real spacecraft angle, z l'
Fig. 4 shows the response of the perfectly matched state synthesiser angle estimate to a step disturbance acceleration, u d , of 80 arcsec
For these tests, an algorithm for calculating the running r .m.s. value of each component of the state estimate is utilized in the simulation program.
s-2, applied at t = 5.52s, the iteration interval having been increased to O.72s. The purpose of this run is to illustrate the capability of the state synthesiser to function correctly with iteration periods which are comparable with or even longer than the periods of natural oscillations of a pla nt with complex poles: impossible to a c hie ve with recursive algorithms due to the stabilit y boundaries imposed on th~ closed loop-poles by the sam pling process. In this case, the period of the flexure mode is about 0.5s, so that there are nearly 1.5 cycles of oscillation between state synthesiser updates.
The first two row s of the table show, respectively, the results with jet noise alone and gyro noise alone. The former c ase shows pessimistic results, since the random plant signal is applied ever y i tera tion but, in practice, there may be many iteratio ns, particularly during limit c ycling periods, in which there are no jet firings. The third row shows the results obtained with simultaneous gyro and jet noises. Subsequent rows of th e table show results obtained with the same levels of gyro and jet noise, but with the introduction of redundant measurement and control samples, the weighting matrix being chosen for least squares fitting, according t o equation (22).
Fig. 4 also sho ws clearly the precise operation of the state synthesiser up to t = 5.52s and the transient erro rs for the 4 iteration s following applic ation of the disturbance ac celera tion step, prec ise estimates being obta ined on the 5th a nd subsequent iterations. The run also serves to illustrate the correct functioning of the state synthesiser with the c ontrol, u, switching before and during the acquisition transient.
As expected, the angl e and less sensitive to jet measurement noise, due integra tions between the jet output.
Sensitivity to plant and measurement noise
rate est imates are noise, than to to the kinem a tic noise and the gyro
The angle estimates contain much less noise than the gyro angle output, since, even for the minimum length synthe siser, five successive measuremnts are used to c onstruct the angle estimate. The minimum length synthesiser thus has inherent filtering properties. The angl e
To provide a measure o f the state synthesiser stochastic performanc e, Table I shows the r.m.s. values of the state estimates obtained with random pla nt a nd measurement inputs
TABLE 1 Sensitivity of State Synthesiser to Plant and Measurement Noise
l\,£ASlREIv£NT
/CINrn:L SEQ.JEN:E LEl\C11-l, s
R.M.S. G'YR:)
R.M.S.
R.M.S. NJISE CN STATE ESTI'v¥\TE
JET
NJI SE
NJISE
arcsec s-I
arcsec
n
zl arcsec
s-2
nz 2 arcsec
nz 4 ar c sec
nz3 arcsec
s-l
nz5 acrsec
s -1
s-2
I
I
5 (Smin)
0
0.28
0.000097
0.00089
0.0022
I 0 .041
0.43 I
!
I
5 (Smin)
0.27
0
0.0023
0.012
0.031
I
i 0.34
5.6 I I
i 5 (Smin)
0.27
0.28
0.0030
0.013
0.033
7
0.27
0.28
0.0021
0.0039
0.013
9
0.27
0.28
0.0018
0.0038
11
0.27
0.28
0.0013
0.0029
i
! II
0.37
1
6. 1
i
0.12
1.1
0.010
0.12
0.80
0.010
0.10
I
._
0.55 '-._
--
200
S. J. Dodds
and rate estimates, in fact, are probably good enough, for many applications, to enable the state synthesiser to be used alone for state estimations. CONCLUSIONS The state synthesiser provides an alternative to a conventional state observer or state estimator and offers the following advantages: (a)
There is no theoretical limit on iteration period.
(b)
In absence of stochastic disturbances and/or plant parameter mismatches, an exact state estimate can be obtained in the minimum possible time, for a given iteration period.
In cases where stochastic performance with the 'minimum sequence' state synthesiser is unsatisfactory, redundant measurement and control samples may be utilized with the introduction of the weighting matrix , W , but with the penalty of increasing the number of arithmetic operations. This is not the case with the recursive estimators in which the gains are adjusted without increasing the complexity of the algorithm. The state synthesiser, however, still produces the theoretically exact state estimate in a finite number of iterations, equal to the number of past measurements used per iteration. In applications where computation speed limitations necessitate iteration frequencies which are comparable with, or even smaller in magnitude than any of the plant poles, a conventional, recursive, estimation algorithm cannot be used and the state synthesiser offers a possible solution. Realistic plant parameter mismatches can be tolerated, but it is recommended that, as a general rule, frequencies of natural vibrational modes of the controlled dynamics are purposely over-estimated in order to avoid potential instabilities. The state synthesiser provides a means of substantially reducing the sensitivity of high precision spacecraft attitude control systems to high level transient disturbance torques within the limitation of the maximum available control torques and the selected control laws. REFERENCES Dodds, S.]. (I981 a). Adaptive, High Precision, Satellite Attitude Control for Microprocessor Implementation. Automatica, Vol. 17, No.4, pp. 563-573. Dodds, S.]. (I981 b). Bang-bang Control Law for Single-Input Time-Invariant Plant. lEE Proc., Vol. 128, Pt D, No.5, pp. 227-232.
State Synthesiser
rr==
201
- --- . - ---- - - ~- · ~ - -- - ~--I
!
I~P~RE~D~'~C~T!~N~G-S~TA~T~E-S~Y~N~TH~E~ . 5~'~SE~R------'
11
I:
[u ]
CC)~TROL VECTOR
I
I
5U
!
~t1ATRIX
ST O .... E
~
X
I!
DETERMINISTiC
DISTl;RBANCE
11
seJ<;OR NOISE
II I 1
I 1
" I,
NOISE
l
CONTROL LAW
SENS ORS
I1
PLANT
DYNAMiCS
'PLANT)
rI
11
:i
---
Jr,
11
k11
Gr..I :-J .".ATRIX
K
IL
OYNAMICS
Fig. I
I
J
11
MODEL
Kt\lMAN
Y
I
==-=X==J
FlU ER
PLANT STATE
Block diagram of control system embodying state synthesiser
--~---- - - -- - -
PLOT
SCA LIN G FACTOR, S
- - ----~----------~---
oS 20
jS
100
arcsec
aroec 5- 2
--~--- -- ~~ --~
u
200
arcsec 5- 2
u
u
~
1- - --- -,
II - T-\.2 - ::
~.~
__
~
__ _ . - - f - ( ! )
'\,-'-----+----+--+-+----- ----+---l-
' - ---l
- r ---
r-~~----"" ~ll..
\\
f--------- -~--
1
~- ---T~-,__ · ~--,_-··,- -________,f---,--+--,----,t-,----
J2.4
). ~ _
__)~~ J
4.'
4.8
5.4
I 7.8
-s
L_::) ____ ::x: ~
8.4
l'lHE (5)
-15
Fig. 4
State synthesis er angle estimate with long iteration period
9.0
9.6
il
202
s.
J. Dodds
(b)
(c)
" -FS
(d)
-"
" -FS
(e)
(g)
-"
(h)
-" -'S,
,
u
"
"
(i) Tt,. I S)
Full-scale (FS) values: 21 = 100 arcsec,
22 = 50 arcsec 5- 1,
Fig. 2
u
200 arcsec s-2,
ud
Control loop responses to step disturbing acc eleration
100 arcsec s-l
S
T1Mf: (S)
0
-2
s- I
5- 1
4200
• 100
- 25
S
'2 5
-]
"
. ]
-7 '1
0
- I
•) 5
l---TI -/(1
n
' ]0
~
.~
1. 2 lid
I 1.8
'"
~(:
en
rt III rt
Id"
(t)
en
'<
2.4
~
rt :T (t)
....'" '"'1 (t)
1.0
1.6
4.1
Fig. 3
Detailed state synthesiser response to step disturbing acceleration IV
o
w
© IFAC Automatic Contro l in Space
~oo rdwijke rhout.
The
~rtherlands
1982
STATE SYNTHESISER; A DIGITAL OBSERVER FOR SPACECRAFT ATTITUDE CONTROL SYSTEMS S.
J. Dodds
Attitude Control Systems Group, Satellite Ditq·sioll , Marconi Space and Defence Systems Ltd., Browns Lane, Th e Airport , Portsmouth , Hampshire , UK
Abstract. The state synthesiser is a form of digital state observer which is optimised for speed of response rather than minimum noise content. The serious limitation often imposed by a Kalman filter on the speed of response of a high preCision control loop to transient disturbance torques (such as produced by orbit change manoeuvres) is removed by switching from the Kalman filter output to the state synthesiser output when the pointing error exceeds a given threshold. In the state synthesiser, the state transition equation for the dynamics is repeatedly applied to construct an estimate of the present state given a finite, minimum length, sequence of past measurements together with the control function operating during the sequence. Simulations are given of a control loop embodying a statE' synthesiser, steady state Kalman filter and sub-time optimal control law (which is adaptive to a disturbing torque estimate) for attitude control of a spacecraft with a flexible structure. A lumped parameter dynamics model is utilized. The principle is not limited to the above application and may be relevant to other time invariant, linear plants of arbitrarily high order. Keywords. Adaptive control; Aerospace computer control; A tti tude control; Bangbang control; Digital computer applications; Discrete systems; Observers; State estimation. perfect performance is obtained in absence of plant and measurement noise or parameter mismatching.
INTRODUCTION The principle function of the state synthesiser is to provide a state estimate which follows, as closely as possible, a rapid Change in the plant (dynamics) state produced by an unexpected disturbance of unknown magnitude. This enables an automatic control loop to respond rapidly, reducing transient errors.
The initial motivation for development of the state synthesiser is the need to minimise transient pointing errors in high precision spacec raft, due to disturbing torques occ urring during firing of orbit adjustment thrusters. Such a problem is encountered on the X-ray astronomy satellite, Exosat. The latter spacecraft employs an alternative scheme embodying gain switching of a more conventional state estimator as described by Dodds (1981 a). Recognition of the dangers of impairing the stability of high order control loops by employing such gain switching, however, leads to the concept of running two state estimators concurrently, one with stochastically optimised gains for fine pointing and the other with higher gains set for fast response. A sudden disturbance causes the optimal state estimator error (the difference between the measurements and corresponding quantity derived from the estimator dynamics model) to grow rapidly at first. The passing of the error beyond a given magnitude threshold is used to indicate the onset of the disturbance and the state input to the control law is switched from the optimal estimator to
The concept is eminently suitable for application where the control law requires an estimate of the complete state vector. Additional state variables which characterize the disturbance are also estimated, as a necessary part of the process of obtaining the plant state estimate and these may also be utilised by the control laws. The controlled plant is assumed time invariant, or, at most slowly time varying but may be of arbitrarily high order, have multiple control inputs and no restrictions on open loop pole locations. The state synthesiser functions in a similar fashion to a finite impulse response filter and employs a non-recursive algorithm in which the estimate of each state variable is obtained as a weighted sum of a finite sequence of past measurements and corresponding control signals. As with the Kalman filter, mathematically
A: ':
62 _ :' "
191
s. J.
192
the 'fas t' e s ti ma to r as lo ng a s the stat e e s tim ato r e rr o r magni t ude rema in s a bove t he thr esho ld, enab li ng con tr o l t o rqu es t o be prod uced rapid ly to combat th e d is t urba nce . Such a sc heme is il lus tr ated in F ig . I. Inc rease o f t he speed of respo nse o f a d ig ita ll y imp leme nted, disc re t e es ti ma to r by inc rease o f ga in s is u lti mately li m ited by t he it e ra t io n pe ri od , h, uns t able osc ill a t io ns o f th e s t a t e est ima t e occu rr ing a t a f r e qu e ncy o f 1/2h Hz if t he ga in lies outside a ce rta in st ab ilit y boundary : This is due to closed loo p po les o f the es t ima to r becoming compa rable in magni t ude wi t h I /h, as th e gains a re inc reased . The sta t e svnthesise r overcomes th is limi t a t ion and also enab les an exact state es t imate to be obtained with a r bitra r y initia l conditions in t he minimum number of algo r ithm it e ratio ns.
Dodds The disturban c e, V( t), is assumed to be deterministi c in na ture and c omposed of c o m ponents whi c h can be ap pro ximate d by so lutions to the st a t e equat io n
W = QW
(3 )
whe re W is a st a t e vec t o r, o f d im e nsio n, w, a ssoc iated with th e d is turb a nce and Q is a fi xe d matrix c har ac teri z ing the di s turban ce process . E xa m ples inc lude st e p di s turbing t o rques and cyc li c d istu rbing t o rqu e s appr ox imate d by a t ru nca t e d Fou r ie r se ri es . The d ist urba nce vec t o r is then ( 4)
V = PW
where P is a fix e d ma tri x o f di me nsio n v x w.
,A, fas t s t a t e observe r , such as th e state svnthes ise r, is necessarily se nsi ti ve to measu remen t and plant noise a nd is t he re fo re int ended to supp lement an op ti mal st a t e estima to r (Kalman fil te rl. In si t ua t io ns where the measu reme nt and plan t noises a re no t too seve re , however , the s ta te synt hesise r may ent ire ly replace the Ka lman fil te r.
St oc has ti c d is t urbances are no t inc lud e d in eq ua ti o ns (1) a nd (2), since t he de rivati on o f th e st a t e s ynthes ise r a lgo r ithm is de t e rmin is ti c in na t ur e and does no t re qu ire kn o wledge o f plan t and me a sur e me nt no ise. The latte r a re inc luded in the simu la ti ons , however, in o rd e r t o assess se ns iti vi t y t o rea li s t ic senso r a nd co ntr o l ac tu a t o r noise .
The simu la t ion results presented in t his pape r a re fo r a tt i tude control about a s ingle axis o f a spacec raft embodying a sing le , domina nt , f lexu re mode. This re lat ive ly s imp le examp le serves to demonst r ate t he ope r at io n o f th e s t a te sy nt hes ise r , but it mus t be rea li sed th a t t he proced ur es given fo r ca lc ul a t io n o f t he s t ate sy nt hesiser matrices a re ge ne r a ll y applicab le and may inc lude simul t aneous t h reeaxis control with linear inte r- ax is co up ling, many f lexu re modes a nd m ulti ple con tr ol ac tua to rs pe r ax is .
Fi na lly , e qu a ti ons (1), (2 ), (3 ) and (4 ) may be comb ined so that , in pa rt it io ned fo r m ,
U
+
(M OJ
Y
+
DU (5 a )
a nd in the sta nd ard fo r m, DEFI NITI O"-! OF PLAN T AND DI STUR BANC E PARAME T ERS In view o f the app licab ili ty t o a ny ti me in va ri a nt, linea r plan t, th e concep t is deve lo ped in ge ne ra l t e r ms , app li ca ti on t o a spec ifi c spacec ra ft dy nam ics being au t oma ti ca ll y included in this general fo r mula t ion. The plant to which the state sy nthesise r app li es !s assumed to be of t he fo r m
Z = FZ
{
x
AX
+ BU
Y
ex
+
(5b) DU
where X is an ove ra ll sta t e vec t o r of d im e ns ion N = n + wa nd A, B a nd C ar e pla nt, inp ut a nd measu reme nt ma tr ic es co rres ponding t o th e pa rtiti oned ma tri ces of e quation (5a). DEVElOPME NT OF THE C O NC EPT
+
GU
+
EV
whe re Z is th e pla nt st ate vec t o r , of di me ns io n n, U is th e co ntr o l vecto r o f d imen sio n r, V is th e di s turban c e vec t o r of d im e ns io n v, F is the plant matri x, B is th e input matrix and E is t he distur banc e ma tr ix. Th e me a sur e ment equ a t io n is of the fo r m Y = MZ + DU
( 2)
where Y is the measure ment vector of dimension m, M is the outp u t matrix and D is the feed-forward matrix.
In o rd e r for the c ur re nt s tate to be estimated in minimum time, witho ut pr ior kno wled ge of the state, the minimum quantit y of input / output obser vation s of the plant must be ut ilized. In the state synthesiser, these observations take the form of a stored past histor y of y( t) and U( t) , sampled at the algorithm iterat ion interval, h. The essential fa c t which determines the minimum number of samples of (Y, U) is that a sufficiently large set of temporariy spaced measurement coordinates serve as well as spacially separated co-ordinates within . Y for the purpose of re-
State Synthesiser constructing the state vector at a given time. In fact, as proven in the algorithm derivation below, the product of the number of samples of Y and the dimension of Y must be at least equal to the order of the system described by equation (5).
193
In general, C is of dimension m x N and is non-invertible, but it is possible to construct an invertible measurement equation by utilizing more than one sample of Y and relating them to X k by repeated application of the state transition equation, as shown below. Suppose that s consecutive input/output samples are utilized, including the current sample (Y k' lI k _I )·
DERIVA TlON OF THE STATE SYNTHESISER ALGORITHM
Then This present derivation assumes constant values of U between algorithm updates. Generalisation to known U(t) between algorithm upda tes, however, should be possible, but is not treat e d in this paper. The state transition equation and measurement equation corresponding to the continuous state equations (5) is Xk '"
(6a)
Y k '" CX k + DU k _ 1
(6b)
Y k-q
=
f
h
o
+
(9)
DlIk_q_1
and inverting the state transition equation: X _ '" ~-I k q
Xk_q+ I _~-I
'V
lIk_q
'" ~-2
X k _q +2 - ~ -2
'f
lIk_q+ I
- ~ -I
'f
(10)
lIk_q
Repeated application of this substitution yields
where cP is the state transition matrix and 'V will be termed the input transition matrix, so that 'I' (h)
CX k _q
)('_q '" "K
J,
'!'
-q X
k
q L: ,;, - i .\'t', I Ik-q + i-I 'I' i =I
-
(11 )
q '" 0, I, 2, .•... S-I
Substituting for Xk q equation (11) Yields -
The state synthesiser algorithm is derived simpl y by repeated application of equation (6) to past samples of (Y, U).
in
equation
(9),
using
q C ~-q ~
Now, it is assumed, as is usually the case, that m < N where m is the dimension of Y and N is the dimension of X. for m '" N, the state synthesiser algorithm is trivial and consists of the inverted measurement equation
Z i=I
+
D lIk _q _ 1 (12)
q '" 0, I, 2, ..... S-I
(8)
(sm x I) ~
Yk Yk - I
Now, equation (12) may be written for every sample, (y k-q' lI k _q _ I ), in partitioned matrix form as follows:
(sm x N) C
----
X k
C~-I
-D
CqJ-I \11
0
0
-D
0
~
.............. , "-
Yk - 2
C~-2
I I
I I I
I t
I
I
I
Yk - s + 1
C~-2
I
I
C~
-(s-I)
(sr x I)
(sm x s r)
~
I
If'
Cq,-I'ji
,
"-
"-
-D
'-"-
"-
cf(s-l)~ Cf(S-2)~
"-
"-
"-
"-
"-
"-
"-
..•..•. '
"-
"-
'--
"-
"-
"-
0
lI _ k 1
0
lIk_2
,
,
0 "-
Cfl~
"-
"-
-D
Uk _ s
(I3)
s. J.
194
The expressions written in parenthesis over equation (13) indicate the matrix dimensions. Now, equations (13) represent s.m linear simultaneous equations in the N unknown components of Xk , and can be solved as they stand only if s.m = N. If s.m < N, there are insufficient equations, tantamount to the (Y) samples containing insufficient information for reconstruction of the state vector. Hence the minimum number of samples can be expressed as
Dodds amount equal to the time taken by the control law to calculate and issue U to the control k actuators. This problem may be avoided by ,A
applying the state transition equation to X in k order to predict the state one iteration step ahead in time. Thus
(18) (14)
smin = int(N/m)
where int(x) is defined as the smallest integer which is ~ x. [n many cases, N/m is not an integer, in which case equations (13) are overdetermined for s = smin. This problem may be overcome, however, by lntroduclng a weighting matrix, W, of dimension N x s.m. First, equations (13) are written in compact form as
I
I
I
issue
of
the
latter
to
the
control
actuators being delayed until sampling instant k+ I. may
be
re-written
in compact
Where [UJ k has been re-defined to include U as its first
r elements.
Sand S
y
u
k are the
fixed, predicting state synthesiser matrices of dimension N x s.m and N x (s+ 1).r, respectively. In terms of the matrices of equation (18),
(Nxsr)
T'I) TI) I
(19)
(15)
where the matrix dimensions are indicated in parentheses. Both sides of equation (15) are then pre-multiplied by W to yield (NxN)
Uk+ I,
~ srxI )
[yJ k = [C]Xk + [DJ[UJk
(Nxsm)
step is allowed for computation of Xk+ I and
Equation (18) form as
( smx s r ) ( smx I ) ( smxN)
Equation (18) constitutes the algorithm of the predicting state synthesiser which is the version recommended for use in automatic control applications. In this case, one iteration
( 1" ' 1 )
(20)
w[yJ k= [W[CJjXk + [WeD] ][U]k
(16)
in which the matrix, W [C] , may be inverted. [t can br shown that for observability of X, Rank [C J = Nand it follows that W may be chosen so that W may be inverted.
Cc J
The basic state synthesiser algorithm is then obtained from equation (16). The synthesised A
state is denoted X in order to distinguish it /\
from X, the symbol usually used to denote the state estimate from a Kalman filter. Then, prfmultiplying both sides of equation (16) by
(21) In absence of plant parameter mismatch and stochastic disturbances, the predicting state synthesiser is capable of observing the state perfectly, after s+ I initial iterations, occupied by an acquisition transient, analogous to the behaviour of a finite impulse response filter. The algorithm is particularly simple in that it is non-recursive and the estimate of each state variable is obtained independently as a weighted sum of a minimum length sequence of past measurements and control values.
LW[CJ] -I yields FIL TERING STATE SYNTHESISER
For
the
cases
in
which
smin
<
N/m
(see
equation (14)), there are m.smin - N redundant ~
Now, Xk is the estimate of the state at the instant of the current measurement sample, y k' and since the next control, Uk, is calculated by the control law as a function of Xk , the control would be applied late by an
measurement coordinates. Suitable choice of W then enables the redundant measurements to be used to advantage, giving the state synthesiser filtering properties. In fact, there is no theoretical upper bound on the number of samples, s, of (Y, U) and in cases where stochastic performance of the minimum
State Synthesiser sequence
(s
smin)
state
synthesiser
is
unsatisfactory, a trade off may be made by increasing s at the expense of increased acquisi tion time «smin+ I)h increased to (s+ I)h). It may further be shown that W may be chosen so that for any s < N / m, the state synthesiser becomes a least squares estimator in the following sense: Let Yk(t) = [C]xk(t) where tk (s+ I)h ~ t ~ tk and Xk(tk) ~ Xk, then if (22)
then Y k(t) is a least squares fit through the measured values Y i' i = k-s, k-s+ I, ..... k-l, k, as
generated
assuming the k values U i , i=k-s+l, k-s+2, ..... k-I, k.
xl x2
from
I:
X
0
0
x3
0
0
0
x4
0
0
_ w
Xs
0
0
0
and
z1
w
control
0
0
0
0
0
2
0 0
0
2
SPACECRAFT ATTITUDE CONTROL APPLICATION For the purpose of demonstrating the state synthesiser, a simple spacecraft example is considered, embodying a single flexure mode. The simulation is confined to gas-jet control, about a single axis, of a spacecraft comprising a main rigid body and a significantly flexible and rigidly attached appendage, such as a solar array, the whole craft exhibiting a single dominant vibrational mode. The flexible mode is represented by a second rigid body with one degree of rotational freedom relative t o the main body, the two bodies being restrained about the control axis by means of a linear torsion spring with zero damping. The gas jets are assumed to apply control torques direc tly to the main rigid body. The relevant plant equation, in the companion form, is
"1 x2
0
x31
0
J
0
w
2
0
l+ A 0
z3
0
-I
0
1I
lo
:l
xl
0 1
x3
x2
i
z4
0
0
0
-I
z sJ
0
0
0
0
y
z 1 + nm
(u + npl
1
,,J
l+A z2
0
+
0
0
19 5
:~
I
I x4
1
l xsJ
(23)
s. J.
196
where z I and z2 are
the
Dodds
pointing error and
pointing error rate, respectively, and z3 and z4 are, respectively, the equivalent flexural-mode displacement and rate. The flexural mode natural frequency is w , the control angular acceleration is u and the angle measurement is y. The moment of inertia ratio between the two bodies of the lumped parameter dynamics model is denoted by A (mode inertia to main spacecraft body inertia). The plant and measurement noises are denoted np and nm respectively.
In cases where limit cycle optimisation requires impulses which are shorter in duration than h, the state synthesiser must be developed further to accommodate U k as a known function of time, defined in the interval (t k , tk+ I). No difficulties are anticipated in this development. Where reaction wheels are used for control, however, the digital control law may calculate control torque demands intended to be held constant for each iteration, rendering the present state synthesiser directly applicable.
The deterministic disturbance is assumed to be a constant disturbing acceleration and is modelled in equations (23) as x5 and z5. Note that the state in X-space is to be estimated, since the control law utilized for the simulation is a bang-bang control law which requires, as input, the state of the plant in the companion form. This control law, which applies to any time invariant, single input, linear plant is presented in detail by Dodds (J98Ib). The state transition equation of the plant is
h
XI
( 1 - c) /w 2
si",
Digital simulation results of a spacecraft attitude control system of similar form to that depicted in Fig. I, using the dynamics equations (23), are presented in Figs. 2, 3 and 4. The computer utilizes a word length of 39 bits. Angular units are given in seconds of arc since applications to spacecraft of high pointing precision are envisaged.
h-slw)/w 2 Qh 2 +(c_1 )/w 2 )/",2
(
( 1 _ c) Iw 2
x2
0
x3
0
0
c
x4
0
0
-ws
c
0
0
0
0
x5
DIGIT AL SIMULA TlON RESULTS
slw
xI
(tIl 2 • (c - I) 1",2) 1",2
( h - si", ) I", 2
x2
(h_s/",)/w 2
(l_c)/w 2
x3
(l_c)/w 2
slw
s
x4
uk
/w
x5
k,1
k
where h is the iteration step length of the onboard digital processing for state estimation and control. .A lso, c = cos w hand s = sin who In this case, N = 5 and m = 1, so that, by equation ([4), smin 5. Now, in general, the
The nominal values of the constant parameters of the simulation are chosen as below to represent one control axis of a three axis stabilised satellite having a flexible appendage with a significant vibrational mode. Inter-axis coupling is not included.
predicting state synthesiser requires smin + 1
Control acceleration from gas jets
samples due to the last term of equation ([2), but, for this example, D = 0 so that only 5 samples of y and n are required for the
= 100 arcsec s-2
~
minimum length sequence to generate X. The predicting state synthesiser matrices, Sy and Sw are therefore, each of dimension 5 x 5 and are generated using state transition matrices similar to that used in equation (24), but with h replaced by -ih, where i = 1, 2, 3 and 4. Note also that, for this example, the control between iterations is assumed constant. A value of 120 ms is chosen for this example, this being the minimum 'on' time of the control jets.
Flexure mode natural frequency, w=12rads- l Lumped parameter model inertia ratio, 3
,,=
State synthesiser iteration interval, h = 0.125 s. The above control produced by a 0.05 operating at a 1 spacecraft moment of m 2•
acceleration would be Nm cold gas thruster m arm with a total inertia of about 400 kg
197
State Synthesiser Control loop responses Fig. 2 shows the control system responses to a disturbing acceleration step of 80 arcsec s-2 (80% of the available control torque) applied at t = 1.2 s and removed at t = 20 s.
simulation to obtain the gyro angle noise. Performance is not drastically reduced by the presence of plant and measurement noise, the peak pointing errors being of the same order as those for the ideal case of Fig. 2b. At t = 1.8 s, the disturbing acceleration estimate, C'd' jumps as the error magnitude threshold switch (rei. Fig.l) changes from the Kalman fi!1er to the state synthesiser. The variations in Lld due
For comparison purposes, Fig. 2a shows the response obtained with a steady state Kalman filter only, in which case the dynamics of the filter dominates the control loop transient, resulting in relatively large pointing errors, the peak error following application of the disturbance acceleration step being about 500 arcsec and that following removal of the disturbance being about -200 arcsec. Note, that in order to show this pointing transient on the same angle error scale as the remaining responses, in which the state synthesiser is brought into play, an unfolded cylindrical plot is utilized. The final portion of the settling transient following the -200 arcsec peak is of limited interest and has been omitted in order to accommodate the plot on the same time scale as the remaining plots.
towards zero, as the state sy nthesiser is again brought into play, following the removal of the high level disturbing acce lera tion. This is followed by a period of noisy disturbing acceleration estimate until about t = 24s, when the Kalman filter is switched in again.
Figs. 2b through 2i show control loop responses in which the error magnitude switch of Fig. I switches to the state synthesiser when the
Figs. 2d to 2i are noise free runs but with various parameter mismatc hes as desc ribed below.
Kalman filter error, magnitude of 2 arcsec.
zI
zI
exceeds
a
In Fig. 2b, the dynamics model parameters of the Kalman filter and state synthesiser are precisely matched to those of the real dynamics and no plant or measurement noise is included. The peak pointing error magnitude during application of the disturbance acceleration is reduced to about 16 arcsec, in contrast to 500 arcsec without the sta~e synthesiser. The peak pointing error following removal of the disturbance acceleration is reduced to about 5 arcsec in magnitude, in contrast to 200 arcsec without the state synthesiser. This run also serves as a standard of comparison for the simulations described below in which realistic imperfections are introduced. Fig. 2c shows the control system response with perfectly matched parameters as in Fig. 2b, but with plant noise and measurement noise. The plant noise is derived from typical cold gas thruster data and is taken as having a constant spectral density of 0.02 (arcsec s-2)2/Hz, being simulated as Gaussian
to the stochastic disturbance are evident in Fig. 2C. At approximately t = 5.5s, the state estimate from the Kalman filter is again used, as shown by the marked reduction in noise on the disturbing acceleration es timate. The end of the Kalman filter transient can just be seen, as the disturbing acce leration estimate monotonically approaches the correc t value, u . d ~ At about 21s, u d steps abruptly down
In Fig. 2d, th e con trol acceleration assumed in the state synthesiser and Kalm a n filter is overestimated by 10%. Again, at about t = 1.8s, the state synthesiser o utput is used and the initial variations of ""ud a bout t he true value, ud, are due to t he finite acquisition transient of the state synthesiser which occu rs after application of the negative contr ol jet. The difference between the control accelera t io ns applied to the re al dynamics and dynam ics model is auto matica ll y treated as a s tep disturbance acceleration of -1 00 - (-110) 10 arcsec s-2 and is added t o the estimate of ud. A
The variations in u d a fter removal of the real disturbance acceleration are due to the alternate pulsing of the over-estim ated control jets. This variation is drastically reduced when the Kalman filter estimate is used after about t = 25s. In Fig. 2e, the control acceleration assumed in the dynamics models is under-estimated by 10%. The additional disturbance acceleration due to application of the negative control jet
noise with an r.m.s. value of 0.28 arcsec s-2 over a bandwidth of 1/2h ~ 4 Hz where h is the iteration interval of 0.12 s. The gyro noise is assumed to be white in rate and, using typical high precIsion rate integrating gyro data, a somewhat pessimistic spectral density
steady state estimate, ud' is, therefore, lower
of 0.017 simulated
being r.m.s.
remarks apply to this case as those given for Fig. 2d.
value of 0.27 arcsec s-I over a bandwidth of 4 Hz. This rate noise is integrated in the
Figs. 2£ and 2g show runs in which the flexure mode natural frequency, w , assumed in the
(arcsec s-I )2/Hz is taken, as Gaussian noise with an
is then -lOO - (-90) = -10 arcsec s-2.
The
~
than ud by 10 arcsec s-2.
Otherwise, similar
s.
198
state synthesiser and Kalman filter is, respectively, over-estimated by 20% and underestimated by 5%. The control loop response is found to be insensitive to over-estimation of w but very sensitive to under-estimation of w . Runs for over-estimation of more than 20% are not presented since errors of this magnitude are unlikely to occur in practice. Increasing the under estimation by more than 5% produces relatively large peak pointing errors (about 70 arcsec at 5%) and can lead to control loop instability.
The oscillations in
~d
occur while the estimate
of the state synthesiser is being used and comparison with the rate trace, z2' shows that they occur at the flexure mode natural frequency. Comparison between F~s. 2f and 2g shows the oscillations of u to be d approximately anti-phase with the oscillations of z2 for the +2096 mismatch but in-phase for the -5% mismatch. The tendancy towards control loop instability with natural frequency under-estimation is probably associated with these phase relationships. Figs. 2h and 2i show runs in which the modelled inertia ratio, A, is, respectively, overestimated and under-estimated by 10%. These results are similar to those for the natural frequency mismatching, in that similar phase relationships occur between 'C'd and z2.
This is
because the natural frequency of the dynamics model is proportional to 1\. Note that although the -10 % m isma tch of Fig. 2i produces an acceptable peak pointing error of about 35 arcsec, there is evidence of control loop instability in the z2 trace where the Kalman filter is being used just before the removal of the disturbing acceleration step at t = 20s. State Synthesiser Acguisition Transient Fig. 3 shows the individual state variable outputs of the perfectly matched state synthesiser during application at t = 1.56s of a high level disturbance acceleration step, u . d The continuous curves are the real state variables and the 'staircase' functions represent the corresponding state synthesiser outputs. During each iteration of 0.12s duration, the predicting state synthesiser computes the encirled points on the 'staircase' traces of Fig. 3. Before application of the step disturbance acceleration, the state synthesiser produces the correct state estimates. Errors in the state estimates are evident for five iterations following application of the disturbance step, in accordance with the theory: The number of steps occupied by the acquisition transient is equal to the product of the dimensions of the measurement vector (I) and the number of state variables (5).
J. Dodds Even during the transient, the estimates tend to follow the true state variables, except for f!I.
ud
in
which
exaggerated.
the
errors
are
the
most
The unusual transient behaviour
A-
of ud' unlike that of any recursive estimator, such as the Kalman filter, is attributed to the oscillation of the elements of the associated row of the Sy matrix. This oscillation of terms is sequence
-I
,~
-2
basically
,
-3
due
•••••
to oscillations in the
which occur due to the
complex eigenvalues of Ijl. COMPARISON WITH RECURSIVE ALGORITHM With an iteration period of 0.12S, it impossible, in the flexible spacecraft to stabilise the state estimator by choice of gains, due to the imaginary
is found example, suitable poles of
the dynamics model at + 12 rad s-I. All attempts result in an oscillation, at 0.24s period of the state estimates about the true values. The problem is removed in the simulations of fig. 2 by iterating the state estimator twice for every iteration of the state synthesiser, with the dynamics model state transition matrix calculated using an iteration period of 0.06s. Comparable results to those of fig. 2b can be obtained, however, using a state estimator with increased but stochastically sub-optimal gains. This is really an unfair comparison, however, since the state estimator has the advantage of operating at twice the iteration frequency of the state synthesiser: Comparison should be made at the same interation frequency, so that in the spacecraft example simulated, the state estimator, even used as an observer with suboptimal gains adjusted for increased speed of response, proves inferior to the state synthesis er; unsatisfactory, in fact, due to the oscillations of the state estimate about the real state, as described above. A recursive algorithm, such as the steady state Kalman filter with fixed gains, involves fewer arithmetic operations per iteration than the state synthesiser but, besides the iteration period limitation referred to above, also requires a theoretically infinite number of iterations to obtain an exact state estimate commencing with arbitrary initial conditions. The state synthesiser, therefore, compare$ favourably for 'fast' state estimates such as required in the application described in this paper.
State Synthesiser Operation of the state synthesiser with extended iteration period
19 9
simulated as for the run of Fig. 2C, but with zero values of demanded control signal and zero real spacecraft angle, z l'
Fig. 4 shows the response of the perfectly matched state synthesiser angle estimate to a step disturbance acceleration, u d , of 80 arcsec
For these tests, an algorithm for calculating the running r .m.s. value of each component of the state estimate is utilized in the simulation program.
s-2, applied at t = 5.52s, the iteration interval having been increased to O.72s. The purpose of this run is to illustrate the capability of the state synthesiser to function correctly with iteration periods which are comparable with or even longer than the periods of natural oscillations of a pla nt with complex poles: impossible to a c hie ve with recursive algorithms due to the stabilit y boundaries imposed on th~ closed loop-poles by the sam pling process. In this case, the period of the flexure mode is about 0.5s, so that there are nearly 1.5 cycles of oscillation between state synthesiser updates.
The first two row s of the table show, respectively, the results with jet noise alone and gyro noise alone. The former c ase shows pessimistic results, since the random plant signal is applied ever y i tera tion but, in practice, there may be many iteratio ns, particularly during limit c ycling periods, in which there are no jet firings. The third row shows the results obtained with simultaneous gyro and jet noises. Subsequent rows of th e table show results obtained with the same levels of gyro and jet noise, but with the introduction of redundant measurement and control samples, the weighting matrix being chosen for least squares fitting, according t o equation (22).
Fig. 4 also sho ws clearly the precise operation of the state synthesiser up to t = 5.52s and the transient erro rs for the 4 iteration s following applic ation of the disturbance ac celera tion step, prec ise estimates being obta ined on the 5th a nd subsequent iterations. The run also serves to illustrate the correct functioning of the state synthesiser with the c ontrol, u, switching before and during the acquisition transient.
As expected, the angl e and less sensitive to jet measurement noise, due integra tions between the jet output.
Sensitivity to plant and measurement noise
rate est imates are noise, than to to the kinem a tic noise and the gyro
The angle estimates contain much less noise than the gyro angle output, since, even for the minimum length synthe siser, five successive measuremnts are used to c onstruct the angle estimate. The minimum length synthesiser thus has inherent filtering properties. The angl e
To provide a measure o f the state synthesiser stochastic performanc e, Table I shows the r.m.s. values of the state estimates obtained with random pla nt a nd measurement inputs
TABLE 1 Sensitivity of State Synthesiser to Plant and Measurement Noise
l\,£ASlREIv£NT
/CINrn:L SEQ.JEN:E LEl\C11-l, s
R.M.S. G'YR:)
R.M.S.
R.M.S. NJISE CN STATE ESTI'v¥\TE
JET
NJI SE
NJISE
arcsec s-I
arcsec
n
zl arcsec
s-2
nz 2 arcsec
nz 4 ar c sec
nz3 arcsec
s-l
nz5 acrsec
s -1
s-2
I
I
5 (Smin)
0
0.28
0.000097
0.00089
0.0022
I 0 .041
0.43 I
!
I
5 (Smin)
0.27
0
0.0023
0.012
0.031
I
i 0.34
5.6 I I
i 5 (Smin)
0.27
0.28
0.0030
0.013
0.033
7
0.27
0.28
0.0021
0.0039
0.013
9
0.27
0.28
0.0018
0.0038
11
0.27
0.28
0.0013
0.0029
i
! II
0.37
1
6. 1
i
0.12
1.1
0.010
0.12
0.80
0.010
0.10
I
._
0.55 '-._
--
200
S. J. Dodds
and rate estimates, in fact, are probably good enough, for many applications, to enable the state synthesiser to be used alone for state estimations. CONCLUSIONS The state synthesiser provides an alternative to a conventional state observer or state estimator and offers the following advantages: (a)
There is no theoretical limit on iteration period.
(b)
In absence of stochastic disturbances and/or plant parameter mismatches, an exact state estimate can be obtained in the minimum possible time, for a given iteration period.
In cases where stochastic performance with the 'minimum sequence' state synthesiser is unsatisfactory, redundant measurement and control samples may be utilized with the introduction of the weighting matrix , W , but with the penalty of increasing the number of arithmetic operations. This is not the case with the recursive estimators in which the gains are adjusted without increasing the complexity of the algorithm. The state synthesiser, however, still produces the theoretically exact state estimate in a finite number of iterations, equal to the number of past measurements used per iteration. In applications where computation speed limitations necessitate iteration frequencies which are comparable with, or even smaller in magnitude than any of the plant poles, a conventional, recursive, estimation algorithm cannot be used and the state synthesiser offers a possible solution. Realistic plant parameter mismatches can be tolerated, but it is recommended that, as a general rule, frequencies of natural vibrational modes of the controlled dynamics are purposely over-estimated in order to avoid potential instabilities. The state synthesiser provides a means of substantially reducing the sensitivity of high precision spacecraft attitude control systems to high level transient disturbance torques within the limitation of the maximum available control torques and the selected control laws. REFERENCES Dodds, S.]. (I981 a). Adaptive, High Precision, Satellite Attitude Control for Microprocessor Implementation. Automatica, Vol. 17, No.4, pp. 563-573. Dodds, S.]. (I981 b). Bang-bang Control Law for Single-Input Time-Invariant Plant. lEE Proc., Vol. 128, Pt D, No.5, pp. 227-232.
State Synthesiser
rr==
201
- --- . - ---- - - ~- · ~ - -- - ~--I
!
I~P~RE~D~'~C~T!~N~G-S~TA~T~E-S~Y~N~TH~E~ . 5~'~SE~R------'
11
I:
[u ]
CC)~TROL VECTOR
I
I
5U
!
~t1ATRIX
ST O .... E
~
X
I!
DETERMINISTiC
DISTl;RBANCE
11
seJ<;OR NOISE
II I 1
I 1
" I,
NOISE
l
CONTROL LAW
SENS ORS
I1
PLANT
DYNAMiCS
'PLANT)
rI
11
:i
---
Jr,
11
k11
Gr..I :-J .".ATRIX
K
IL
OYNAMICS
Fig. I
I
J
11
MODEL
Kt\lMAN
Y
I
==-=X==J
FlU ER
PLANT STATE
Block diagram of control system embodying state synthesiser
--~---- - - -- - -
PLOT
SCA LIN G FACTOR, S
- - ----~----------~---
oS 20
jS
100
arcsec
aroec 5- 2
--~--- -- ~~ --~
u
200
arcsec 5- 2
u
u
~
1- - --- -,
II - T-\.2 - ::
~.~
__
~
__ _ . - - f - ( ! )
'\,-'-----+----+--+-+----- ----+---l-
' - ---l
- r ---
r-~~----"" ~ll..
\\
f--------- -~--
1
~- ---T~-,__ · ~--,_-··,- -________,f---,--+--,----,t-,----
J2.4
). ~ _
__)~~ J
4.'
4.8
5.4
I 7.8
-s
L_::) ____ ::x: ~
8.4
l'lHE (5)
-15
Fig. 4
State synthesis er angle estimate with long iteration period
9.0
9.6
il
202
s.
J. Dodds
(b)
(c)
" -FS
(d)
-"
" -FS
(e)
(g)
-"
(h)
-" -'S,
,
u
"
"
(i) Tt,. I S)
Full-scale (FS) values: 21 = 100 arcsec,
22 = 50 arcsec 5- 1,
Fig. 2
u
200 arcsec s-2,
ud
Control loop responses to step disturbing acc eleration
100 arcsec s-l
S
T1Mf: (S)
0
-2
s- I
5- 1
4200
• 100
- 25
S
'2 5
-]
"
. ]
-7 '1
0
- I
•) 5
l---TI -/(1
n
' ]0
~
.~
1. 2 lid
I 1.8
'"
~(:
en
rt III rt
Id"
(t)
en
'<
2.4
~
rt :T (t)
....'" '"'1 (t)
1.0
1.6
4.1
Fig. 3
Detailed state synthesiser response to step disturbing acceleration IV
o
w