Äîêóìåíò âçÿò èç êýøà ïîèñêîâîé ìàøèíû. Àäðåñ îðèãèíàëüíîãî äîêóìåíòà : http://num-anal.srcc.msu.ru/list_wrk/ps1/ch2st2.ps Äàòà èçìåíåíèÿ: Tue Dec 17 13:00:30 2002 Äàòà èíäåêñèðîâàíèÿ: Mon Oct 1 20:45:59 2012 Êîäèðîâêà:

Mathematical Modeling of Complex Information Processing Systems 61
MAXIMIN TESTING OF SATELLITE STABILIZATION
V. A. Sadovnichii 1 , V. V. Alexandrov 2 , S. S. Lemak 2 , and D. Vera Mendoza 3
A system for accuracy testing of algorithms now in use for the active magnetic stabilization of
satellite attitude is considered. The problem on accuracy testing of stabilization of equilibrium
positions for a control stochastic system is formulated for the case when some parametric and time­
varying perturbations are known. Some conditions for finding the optimal counterstrategy and the
best result of testing are given. An example of quality testing of the system for the active magnetic
stabilization of pitching oscillations of an artificial Earth satellite is examined.
1. Introduction. At present, orientation systems are widely used to control the attitude of small satellites.
Measuring sensors and operating actuators of these systems interact with the Earth's magnetic field. Such
systems offer certain advantages: high reliability, small weight, ease of manufacture, and relatively low cost
(this is very important for small satellites).
The development of an orientation system includes testing its quality in addition to the synthesis of control
algorithms. The use of computerized land­based testing devices is one of the possible methods for solving this
problem. These devices significantly simplify the processes of production and testing of satellites. Figure 1
illustrates the functional block diagram of one of the possible systems for testing the algorithms of magnetic
stabilization.
Fig. 1
The algorithms intended to control the satellite attitude are represented in Figure 1 as Block 1 to be tested.
The testing system regards these algorithms as a ``black box'' with known inputs and outputs. This system
forms outer and inner conditions under which the tested block operates. Block 2 is a basic one; it is supposed
to model the dynamics of the body. This block performs the integration of equations for the center­of­mass
matrix of the satellite as well as the integration of equations for the motion about this center under the action
of control and disturbance torques. It is possible to include some additional elements into the model: booms for
gravitational stabilization, panels of solar batteries, etc. The satellite is controlled by magnetic coils. Block 1
1 Moscow State University, Moscow, 119899, Russian Federation, e­mail: vsadovn@rector.msu.ru
2 Faculty of Mechanics and Mathematics, Moscow State University, Moscow, 119899, Russian Federation,
e­mail: valex@moids.math.msu.ru; lemak@moids.math.msu.ru
3 Facultad de Ciencias Fisico Matematicas, Universidad Autonoma de Puebla, Apdo. Postal 1364, Puebla,
M'exico, e­mail: vmendoza@fcfm.buap.mx

62 Mathematical Modeling of Complex Information Processing Systems
generates controlling signals for the block of forming the controlling magnetic torques (Block 3). In addition,
Block 3 uses the information on strength of the Earth's magnetic field whose model is realized in Block 5. This
block may have several levels varying in complexity of description of the Earth's magnetic field. The disturbance
magnetic torque acting on the satellite is formed in Block 6. The gravitational forces and torques acting on
the satellite are modeled in Blocks 7 and 8. The gravity field is assumed central when the gravitational torques
are calculated. Some models of various levels of complexity for the description of the Earth's gravitational
potential are used to model the gravitational forces (Block 7). Block 9 is intended to model the behavior of the
atmosphere. Blocks 5, 7, 8, and 9 altogether form the outer conditions that influence the motion of the body
under control.
The inner conditions affecting the body are formed by measurement errors, errors of operating actuators,
intrinsic magnetic torques as well as by a number of malfunctions and emergency situations arising during
the process of operation of the satellite. It is assumed that the tested control system may use the readings
of magnetic and optic sensors. The errors of the sensors are modeled in such a way that their additive,
multiplicative, and random components are taken into account. The statistical description of purely random
components of measurement errors as well as the description of the functional set for nonrandom components
are used in the following two cases: when operation of the sensors is modeled in Blocks 10 and 11 and when
Block 12 (controlling the blocks listed above and realizing the strategy of testing) is active during the process
of determination of the worst perturbations. During the test, Block 12 can change the value of any parameter
influencing the motion of the system, interrupt the analysis of the motion of the satellite and restart it with
new initial data, etc. This block has the following input parameters: descriptions of sets of permissible ranges
in specification of the designed orbit as well as descriptions of ranges of the mass and geometric parameters
of the satellite and its elements (booms, solar batteries, etc.). Descriptions of sets defining the structure and
parameters of each component of measurement errors of the available sensors are also specified. In addition, a
description of atmospheric changes and data on gravitational and magnetic perturbations are delivered at the
input of Block 12.
2. Formulation of the problem of testing. The testing procedure realized in Block 12 is based
on obtaining the reliable data that characterize the stabilization accuracy for the case when the parametric
and time­varying perturbations preventing the system from being stabilized are the worst possible. Such an
approach was proposed in [1], where formulations of maximin problems provide the basis for the testing strategy.
In this paper, this approach is extended to the case when a statistical description is known for a certain part of
perturbations, which may be represented, therefore, by steady stochastic processes with a given spectral density.
As is known, such perturbations can be modeled by perturbations at the output of a linear stochastic system
whose input actions may be regarded as white noise with a given intensity [2]; without loss of generality, hence,
we may assume that the random perturbations acting on the system are modeled by the vector white noise with
a known intensity. As shown below, in this case the testing procedure [1] leads to game problems for dispersion
relations of a stochastic system (these relations are expressed in terms of deviations from the required motion
of the system).
Suppose the equations
—
y = ~
f(y; w; v; ¸)
describe the motion of a control system. Here w is the vector of control signals: w(\Delta) 2 ~
W = fw(\Delta) 2 L s
2 j w(t) 2
\Omega ae R s g; v is the vector of perturbations influencing the stabilization accuracy adversely: v(t) 2 V = fv(\Delta) 2
L l
2 j v(t) 2 V ae R l g; ¸ is the vector of random perturbations modeled by a steady stochastic process.
Let us consider the system of equations
—
x = f(x; v; u; ¸) (1)
given in deviations x = y \Gamma ~
y from the desired motion ~
y(t) belonging to the preassigned trajectory set Y =
f~y(t) 2 Y 0 ae R n ; t 0 6 t 6 t k g. Here u = w \Gamma ~
w are the unutilized resources in control for solving the
stabilization problem. Usually, the control signals are formed as functions of the time t and of the deviation
estimates ~ x: u = u(t; ~
x). The accuracy of robust stabilization is defined by the functional class
J = M
2
4
t k
Z
t 0
(x ? Gx+ u ? Nu) dt + x ? (t k )Sx(t k )
3
5 (2)
where the matrices G and S are constant, symmetric, and positive semidefinite (G ? = G ? 0, S ? = S ? 0)
and the matrix N is positive definite (N ? = N ? 0). It is assumed that t 0 is a fixed instant of time and t k is
the instant when the given manifold is first reached (in particular, t k may be fixed).

Mathematical Modeling of Complex Information Processing Systems 63
The best accuracy J 0
0 (or the best result in terms of accuracy of robust stabilization) is defined as
J 0
0 = sup
v2V
inf
u(\Delta)2U
J (3)
The best accuracy concept coincides with that of minimax stabilization if the upper and lower bounds are
attainable and a saddle point (v 0 , u 0 ) is in existence:
J 0
0 (u 0 ; v 0 ) = max
v2V
min
u2U
J(u; v) = min
u2U
max
v2V
J(u; v) (4)
The vector function v depends on the deviations and control values: v = v(t; x; u). Thus, there exists a certain
discrimination of control; in accordance with the Krasovskii theorem [3], differential game (3) has the saddle
point (4) (the instant t k is fixed).
The testing procedure may be decomposed into the following three stages.
The first stage (preliminary). The best accuracy J 0
0 (a value of game (4)) and the optimal strategy v 0 (\Delta) of
testing are determined.
The second stage (basic). The process of testing is simulated with the known strategy v 0 (\Delta).
The third stage (final). Two comparisons are performed: a) the best accuracy J 0
0 with the point and (or)
interval estimates of the functional J(u r ; v 0 ) obtained in accordance with the results of simulation; b) the
stabilization strategy u r (\Delta) with the best control u 0 (\Delta) obtained at the first stage.
3. Determination of the optimal counterstrategy and a value of the game. Let us consider
the first stage of the testing procedure: obtaining the best accuracy and forming the optimal strategy. To
accomplish this, as was stated above, maximin problem (4) should be solved.
The following system of equations in deviations at the first approximation is introduced:
—
x = A(t; p) x + B(t; p) u +C(t; p) q +D(t; p) ¸ (5)
Here ¸ is the vector white noise with a known intensity matrix F (t) ? 0 (M [¸(t)] j 0, M [¸(t) ¸ ? (t 0 )] =
F (t) ffi (t \Gamma t 0 ), F (t) is a positive semidefinite symmetric matrix, and M is the expectation operator); M [x(0)] = e 0
is the mean value of the deviations x(0) and M
\Theta ffi
x (0) ffi
x ?
(0)
\Lambda = P 0 is the covariance matrix; p and q are the
components of the vector v = fp; qg of perturbations being varied. The subvector p describes the parametric
perturbations acting on the control system. It is assumed that the perturbations p are constant in time and
belong to a known set p 2 W ae R k . The perturbations q are regarded as arbitrary piecewise continuous
functions whose values belong to a given set q(t) 2 Q. In the problem of magnetic stabilization, for example, the
parametric perturbations p are errors of initial orientation of the instrumental frame and the perturbations q(t)
are the ``parasitic torques'' arising because of intrinsic magnetization of the satellite.
Equation (5) should be understood as a brief formal representation of the corresponding integral equation,
whereas ¸(t) should be understood as the generalized derivative of the normal process i(t) with independent
increments (¸ dt = di). The control u belongs to the set “
U of stochastic processes that admit representation (5).
In particular, this representation is possible [4] if u = K(t)x, where t 2 [t 0 ; t k ], K(t) 2 K, and K is the set of
functions whose values are (n \Theta s)­matrices with piecewise continuous elements.
The problem
J 0
0 = max
p2W; q2Q
min
u2 ~
U
J(p; q; u) (6)
is further solved for auxiliary system (5) and the extended set u(\Delta) 2 ~
U (U ae ~
U ae “
U ) of controls. Here
~
U = fu(\Delta) : u = K~xg and K(t) is a matrix with bounded piecewise continuous elements. During the process
of solving this problem, the optimal counterstrategy p 0 , q 0 and the auxiliary stabilization algorithm u 0 are
determined. The algorithm u 0 , which provides the best accuracy of stabilization, is found from the inner
minimization problem.
If complete and accurate data on the deviations x are available and the control resources are unlimited
(\Omega = R s ), then the auxiliary optimal control u 0 has the form of linear feedback in deviations: u 0 = K(t)x.
If only incomplete data are available, then the auxiliary stabilization algorithm is formed with the use of
the deviation estimates ~
x produced by the estimation system on the basis of data from measuring sensors. These
data are represented by the measurement vector z = Hx+H 1 q + j, where j are the instrument errors of the
measuring sensors being described by the white­noise process [2]. It is assumed that the estimation algorithm
has the form [2]
—
~
x = A(t; p) ~
x +B(t; p) u + C(t; p) q + ~
K(z \Gamma H~x)
z = Hx+H 1 q + j
(7)

64 Mathematical Modeling of Complex Information Processing Systems
where ~
K(t) is the feedback matrix in this algorithm.
The stabilizing control is determined in the form u = K~x, where ~
x are the deviation estimates obtained by
estimation algorithm (7).
By the separation theorem [4], the problem of finding the best accuracy (6) for system (5), (7) may be
formulated as follows:
J 0
0 = max p;q
J 0 (p; q); J 0 (p; q) = min
K; ~
K
M
2
4 x(t k ) ? Sx(t k ) +
t k
Z
0
(x ? Gx + u ? Nu) dt
3
5
Here K and ~
K are the feedback matrices in the problems of optimal stabilization and estimation, respectively.
If q j 0, then the inner minimization problem is reduced to that of choosing a linear stochastic controller
on the basis of the quadratic criterion with incomplete measurements.
Let us consider the general case q 6j 0 and assume that complete data on the deviations x are available. If
the control is given in the form of the linear feedback u = Kx, then performance criterion (2) can be rewritten
as
J = Tr
2
4 SP x (t k ) +
t k
Z
0
(GP x + K ? NKP x ) dt
3
5 (8)
where P x is the matrix of second torques (P x = M [xx ? ]) and Tr is the operation of calculating the matrix
trace: Trka ij k =
n
P
i=1
a ii (the sum of the diagonal elements).
The matrix P x satisfies the dispersion relation
—
P x = (A +BK)P x + P x (A + BK) ? + Cqe ?
x + e x (Cq) ? +D ? FD; P x (0) = P 0
x (9)
The stochastic problem of finding a minimum of functional (2) is reduced to the deterministic extremal
problem
J 0 = min
K
J
under condition (9), where the elements of the matrices K and P x are regarded as controls and phase coordinates
of dimension n(n + 1)
2 , respectively.
Note that P x = P + e x e ?
x (where P = M
Ÿ
ffi
x ffi
x ?
–
, ffi
x= x \Gamma e x , e x = M [x(t)]) and e x satisfies the equation
—
e x = (A +BK) e x + Cq (10)
Then it follows from expression (8) that
J 0 = Tr
2
4 SP x (t k ) +
t k
Z
0
(GP x + K ? NKP x ) dt
3
5 = J 0
1 + J 0
2
J 0
1 = Tr
2
4 SP (t k ) +
t k
Z
0
(GP +K ? NKP ) dt
3
5
J 0
2 = e ?
x Se x +
t k
Z
0
e ?
x (G +K ? NK) e x dt
(11)
The solution of extremal problem (6) J 0
0 = J(p 0 ; q 0 ) is obtained by a method of successive approximations.
The inner minimization problem with respect to the coefficients K should be solved at each step of the
iteration procedure. We solve this problem by the method of successive approximations, too. The feedback
coefficients are represented by two components as follows: K = K x + K e . The first component K x (K e should
be fixed) is obtained as a solution of the problem min
Kx
J 0
1 ; therefore, K x = \GammaN \Gamma1 B ? L, where L is the symmetric
matrix satisfying the differential equation
—
L + LA 0 +A 0? L + G 0 \Gamma LBN \Gamma1 B ? L = 0; L(t k ) = S (12)

Mathematical Modeling of Complex Information Processing Systems 65
and A 0 = A +BK e , G 0 = G+K e NK ?
e . The funcional J 0
2 can be represented as
J 0
2 = e ?
x Se x +
t k
Z
0
(e ?
x G 00 e x u ?
e Nu e ) dt
where u e = K e e x and G 00 = G + K ?
x NK x . The coefficients K e are obtained as a solution of the minimization
problem for the quadratic funtcional J 0
2 and linear system (10):
J 0
2 = e ?
x Se x +
t k
Z
0
(e ?
x G 00 e x u ?
e Nu e ) dt ! min
ue
—
e x = A 00 e x + Bu e + Cq
(13)
where u e = K e e x , G 00 = G+K ?
x NK x , and A 00 = A + BK x . After several iterations, we obtain optimal values
of the coefficients K as a solution of the inner problem.
Now we turn to the outer problem. If p 2 W is fixed, then the problem of finding max
q(\Delta)2Q
J 0 is reduced to
the extremal problem of control in q(t) in the linear system
—
e x = A(p)e x + B(p)Ke ~ x + Cq(t) (14)
with the quadratic performance criterion
J 0
2 = e ?
x Se x +
t k
Z
0
e ?
x (G + K ? NK) e x dt ! max
q(\Delta)2Q
(15)
Its solution may be obtained, for example, by the Krylov--Chernous'ko method of successive approximations [5]
whose convergence is proved in some special cases (see [6]). At the next step (once the perturbations q(t) have
been determined), the problem of finding max p2W
J 0 may be considered as a problem of nonlinear programming in
parameters p 2 W , which can be solved by an appropriate numerical method [7].
4. Accuracy testing of the stabilization algorithm. At the second step of the strategy, the stabilization
algorithm of forming control signals u r is tested with the use of the results obtained at the first step. This
algorithm may be unknown for the tested system, i.e., the system considers this algorithm as a ``black box''
with known inputs and outputs.
The testing procedure requires statistical tests for system (1) to be performed and estimates the performance
criterion (2) with the use of methods of mathematical statistics. For this purpose, system (1) is integrated with
known values of p 0 and q 0 (t); to accomplish this, the errors ¸ in the system and the readings of the sensors
are modeled. The control signals u r are obtained at the output of the ``black box'' and system (1) is closed.
After processing the results of integration, the statistical mean value of stabilization accuracy of the algorithm
J 0 (p 0 ; q 0 (t); u r ) is calculated; this value is compared with the best accuracy J 0
0 obtained at the first stage:
Ÿ = J 0 (u r )
J 0
0
This comparison allows us to judge the quality of stabilization by the control signals u r .
The interval estimates of the functional J 0 (p 0 ; q 0 (t); u r ) permit determining the number of necessary sta­
tistical tests at the second stage. If the best accuracy J 0
0 falls on the confidence interval in J 0 (p 0 ; q 0 (t); u r ),
then we should return to the second stage and increase the number of tests.
5. Testing of algorithms for active magnetic stabilization. 1. Let us consider the application of the
above algorithm to the problem of testing the algorithms for active magnetic stabilization of the satellite. The
equations
A¨ff 1 + ! 0 (C +B \Gamma A) —
ff 3 + 4! 2
0 (B \Gamma C) ff 1 =M a
1 + M b
1 + \DeltaM g
1
B¨ff 2 + 3! 2
0 (A \Gamma C) ff 2 =M a
2 + M b
2 + \DeltaM g
2
C ¨
ff 3 + ! 0 (B \Gamma A \Gamma C) —
ff 1 + ! 2
0 (B \Gamma A) ff 3 =M a
3 + M b
3 + \DeltaM g
3
(16)
describe small oscillations about the center of mass in projections on the axes of the coupled system. Here A,
B, and C are the torques of inertia; ! 0 (t) is the angular velocity of orbital motion; ff = (ff 1 ; ff 2 ; ff 3 ) is the vector
of small rotation of the coupled system Mx 1 x 2 x 3 relative to the orbital system My 1 y 2 y 3 .

66 Mathematical Modeling of Complex Information Processing Systems
Let us introduce the following notations: ffi = fffi ij g is the matrix of rotation of the axes Mz 1 z 2 z 3 of the
instrumental frame relative to the coupled system Mx; u = (u 1 ; u 2 ; u 3 ) are the projections of the active torque
on the axes of the system Mz (this torque is generated by current coils); (I n
1 ; I n
2 ; I n
3 ) are the projections of the
intrinsic magnetic torque arising because of magnetization of the satellite; T = (T 1 ; T 2 ; T 3 ) are the projections
of the Earth's magnetic strength vector on the axes of the orbital system My; \DeltaM g
1 , \DeltaM g
2 , and \DeltaM g
3 are the
projections of the disturbance torques (for which a statistical description is available only) on the axes of the
coupled system (for example, these torques may be represented by the truncated Fourier series
\DeltaM g
i = d 1i + d 2i cos ! 0 t + d 3i sin ! 0 t + d 4i cos 2! 0 t + d 5i sin 2! 0 t + d 6i ; i = 1; 2; 3 (17)
where d 1i , d 2i , d 3i , d 4i , d 15 = const and d 6i is the high­frequency component of the white­noise type).
The projections of the active magnetic torque on the axes of the coupled system are denoted by u x =
[E \Gamma “ ffi ] u z . The projections of the control torque M a = u \Theta T are written down in the form
M a
1 = u 2 (T 3 \Gamma T 2 ffi 1 ) + u 1 (T 2 ffi 2 + T 3 ffi 3 ) \Gamma u 3 (T 2 + T 3 ffi 1 ) +
X
i;j
ff i u j – ij1 (ffi)
M a
2 = u 3 (T 1 \Gamma T 3 ffi 2 ) + u 2 (T 1 ffi 1 + T 3 ffi 3 ) \Gamma u 1 (T 3 + T 1 ffi 2 ) +
X
i;j
ff i u j – ij2 (ffi)
M a
3 = u 2 (T 3 \Gamma T \Gamma 2ffi 1 ) + u 1 (T 2 ffi 2 + T 3 ffi 3 ) \Gamma u 3 (T 2 + T 3 ffi 1 ) +
X
i;j
ff i u j – ij3 (ffi)
(18)
Here – ijk (ffi) are the coefficients at the nonlinear components of the active control torque (18). Similar expressions
can be written for the components of the disturbance magnetic torque.
The measurement vector may be represented in the form z = T m
z \Gamma T \Lambda , where T m
z is the vector composed
of the measured components of the field strength at the point M and T \Lambda is the information on projections of
the Earth's magnetic field obtained aboard the satellite. The equality T \Lambda = T y is valid when the center­of­mass
position M of the satellite is known with certainty. The representation
T m
z = [E + \Delta]
h
E + “ ffi 0
i h
E + “
ffi
i
[E + “
ff] T y + ¸ 0 + ¸ 1 sin ! 0 t + ¸ 2 cos ! 0 t + ¸ 3 (19)
holds up to small quantities of second order. Here ¸ 0 = (¸ 0
1 ; ¸ 0
2 ; ¸ 0
3 ) is the systematic error of measurement
(zero drift of the sensors); ¸ 1 and ¸ 2 are the time­dependent components of the systematic error; ¸ 3 are the
high­frequency error components represented by white noise; ffi 0 is the vector of small rotation of the sensors'
sensitivity axes relative to the system Mz; \Delta is the matrix of relative measurement errors. The errors in the
measurements at the scale of magnetometers are represented by the diagonal elements (\Delta 11 ; \Delta 22 ; \Delta 33 ), whereas
the off­diagonal elements (\Delta ij ) are the coefficients of cross sensibility of the magnetic sensor (these coefficients
arise because of nonorthogonality of the measurement axes along which the magnetic sensors are located).
If we use the projections on the axes of the system Mx, then relations (19) are written down in the form
z k = \Delta kk T 1 + (ff \Theta T ) k +
X
i
ffi i Ÿ 2
ik (T y ) +
X
i
ff i Ÿ 3
ik (T y ) +
X
i;j
ff i T j Ÿ 4
ijk (\Delta ij ; ffi; ffi 0 )
+ ¸ 0
k + ¸ 1
k sin ! 0 t + ¸ 2
k cos ! 0 t + ¸ 3
k
(20)
The additional equations
( —
d 2i = d 3i ;
—
d 3i = ! 2
0 d 2i ;
( —
d 4i = d 5i ;
—
d 5i = 4! 2
0 d 4i ;
( —
¸ 1
i = ¸ 2
i
—
¸ 2
i = ! 2
0 ¸ 1
i
—
d 1i = 0; —
¸ 0
i = 0; —
Ÿ i = 0; i = 1; 2; 3
used to describe the constants in expressions (17) and in the measurement vector (20) are supplemented to
equations (18). Assuming that ff i and ffi j are small, we ignore the nonlinear terms in (18) and (20). Then, the
equations of motion of the control system may be represented in the form
—
x = A(t; p) x +B(t; p) u + C(t; p) q + D(t) ¸
—
~
x = A(t; p) ~
x +B(t; p) u + C(t; p) q + ~
K(z \Gamma H~x)
z = H(t; p) x + H 1 (t) q + j

Mathematical Modeling of Complex Information Processing Systems 67
Here x = (ff 1 ; ff 2 ; ff 3 ; —
ff 1 ; —
ff 2 ; —
ff 3 ; d 1i ; d 2i ; d 3i ; d 4i ; d 5i ), i = 1; 2; 3; u = (u 1 ; u 2 ; u 3 ); p = (ffi 1 ; ffi 2 ; ffi 3 ); q = (I n
1 ; I n
2 ; I n
3 );
¸ = (d 61 ; d 62 ; d 63 ; ¸ 3
1 ; ¸ 3
2 ; ¸ 3
3 ).
This system is similar to equations (5) and (7); therefore, the strategy of testing (i.e., obtaining the best
accuracy and forming the optimal counterstrategy (p 0 ; q 0 )) can be found solving the extremal problem (6).
2. Let us consider in more detail the calculation of the following three characteristics: a value of the
game, the optimal counterstrategy, and an estimate of stabilization quality. For doing this, we consider the
stabilization of pitching oscillations of the satellite moving along the circular polar orbit.
It is assumed for simplicity that the Earth's magnetic field strength is formed by the field T specified by
the dipole and the fluctuating component \DeltaT . Let the satellite move uniformly along the circular orbit. The
orbital system of coordinates is denoted by Oy 1 y 2 y 3 . The magnetic torque M a generated by the current coil is
represented in the form M a = M 0
u \Lambda u, where u = f1; 0; \Gamma1g is the control.
The coordinate system associated with the principle axes of inertia of the satellite is denoted by Ox. Let ff
be the angle of rotation of the system Ox relative to the orbital system Oy. The system Oz is associated with
the instrumental frame of the satellite. This system is rotated relative to Ox through the angle ffi 0 + ffi , where ffi 0
is the required value for the angle of rotation and ffi is the error of realization of the instrumental frame. The
equation for pitching oscillations may be written down in the form
B¨ff + 3! 2
0 (A \Gamma C) sin ff cos ff = M a +M b + \DeltaM ¸
Here A, B, and C are the inertia torques of the satellite; M a and M b are the control and disturbance magnetic
torques; \DeltaM ¸ are the disturbance torques of external forces. Let \DeltaT y = ( ~
\Delta 1 ; 0; ~
\Delta 3 ) be the projections of
perturbations of the Earth's magnetic field strength on the axes of the system Oy. The intrinsic magnetization
of the satellite causes the magnetic torque whose projections on the axes of the coupled system Ox are I oe =
(n 1 ; 0; n 3 ). It is assumed that a single current coil (whose magnetic torque is oriented along the z 3 ­axis) is
mounted on the satellite. Then, the control torque is
M a = 0:5M 0
u u
h
2e 0
~
\Delta 1 cos(ff + ffi 0 + ffi ) \Gamma 2e 0
~
\Delta 3 sin(ff + ffi 0 + ffi )
\Gamma B 0
\Gamma
3 sin(ff + ffi 0 + ffi \Gamma Ü ) + sin(ff + ffi 0 + ffi + Ü )
\Delta\Lambda
where B 0 = ¯ 0 H 0 =r 3
0 , H 0 is the magnetic torque of the Earth's dipole, ¯ 0 is the coefficient of magnetic
permittivity, and r 0 is the orbit radius of the satellite.
The disturbance torque is
M b = B 0 M 0
u
\Theta oe 3 (\Delta 1 + sin Ü ) + oe 1 (\Gamma\Delta 3 \Gamma 2 cos Ü ) + ff
\Gamma oe 1 (\Gamma\Delta 1 \Gamma sin Ü ) + oe 3 (\Gamma\Delta 3 \Gamma 2 cos Ü )
\Delta\Lambda
If we retain only the terms linear in ff, then the oscillatory equations can be written in the dimensionless form
as follows:
¨
ff +
\Gamma b 2 \Gamma b 50 oe 1 \Gamma b 60 oe 3 + b 20 u
\Delta ff = u[b 00 + b 21 ffi ] + oe 1 b 31 + oe 3 b 41 + ¸ 1
+ \Delta 1 [u(b 70 u + b 71 (ffi + ff)) + b 40 oe 3 ] + \Delta 3 [u(b 80 u + b 81 (ffi + ff)) + b 80 oe 1 ]
(21)
Here Ü = ! 0 \Lambda t, M 0 = M 0
u B 0
! 2
0 B , \Delta i = e 0
~
\Delta i
B 0
, oe i = n i
M 0
u
are the dimensionless variables, and
b 2 = 3 A \Gamma C
B
; ¸ 1 = \DeltaM ¸
B! 2
0
; b 10 = 1; b 30 = \GammaM 0 ; b 40 = M 0
b 00 = M 0
`
\Gamma 3
2 sin(ffi 0 \Gamma Ü ) + 1
2 sin(ffi 0 + Ü )
'
; b 20 = M 0
` 3
2 cos(ffi 0 \Gamma Ü ) + 1
2 cos(ffi 0 + Ü )
'
b 30 = \Gamma1; b 40 = 1; b 21 = \Gammab 20
b 50 = M 0 sin Ü; b 60 = 2M 0 cos Ü; b 70 = M 0 cos ffi 0 ; b 80 = \GammaM 0 sin ffi 0
b 31 = \Gamma2M 0 cos Ü; b 41 = M 0 sin Ü; b 71 = \GammaM 0 sin ffi 0 ; b 81 = \GammaM 0 cos ffi 0
The Mexican satellite ``Satex'' travelling along the orbit of radius r 0 = 800 km has the active magnetic torque
M 0
u = 4 A \Delta m 2 ; for this case, the parameters M 0 and b 2 are equal to 1:519 and 2:936, respectively.
Let p = ffi , q 1 = oe 1 , q 2 = oe 3 , ¸ 1 = \DeltaM ¸
B! 2
0
, ¸ 2 = \Delta 1 , ¸ 3 = \Delta 3 , x = (ff; —
ff), q = (q 1 ; q 2 ), ¸ = (¸ 1 ; ¸ 2 ; ¸ 3 ; ¸ 4 ).
Equation (21) can be rewritten up to ff\Delta i , ffi \Delta i , ffoe i \Delta i as follows:
—
x = A(t)x +B(t; p)u +C(t)q + D(t)¸ (22)

68 Mathematical Modeling of Complex Information Processing Systems
Here
A =
` 1 1
\Gammab 2 0
'
; B(p) =
` 0
b 00 + b 20 p
'
; C =
` 0 0
b 31 b 41
'
D =
` 0 0 0 0
b 10 b 20 b 70 b 80
'
; M
\Theta
¸ 4 (t) ¸ 4 (t 0 )
\Lambda
= D 0 e \Gammafl t ffi (t \Gamma t 0 )
The parameters D 0 and fl are chosen in such a way that the solution of auxiliary system (22) on the segment
[0; t k ] differs little from the solution of nonlinear system (17).
It is assumed that the error p 2 W of initial orientation of the instrumental frame belongs to a known
segment and the components of the vector q(t) are arbitrary time functions whose values are known up to the
set q(t) 2 Q. The components of the vector ¸ are regarded as independent white noise with intensities ~
F i .
At the first stage of testing, problem (6) of finding the best accuracy is solved by a successive approximation
method. It is assumed for simplicity that complete data on the state vector x are known.
The worst perturbations q 0 (t) are obtained by solving the linear problem (14) with the quadratic perfor­
mance criterion (15).
The dependence of the performance criterion J 0 (p; q 0 ) on p is presented in Figure 2 for the following initial
data:
t k = 12; ~
F(1; 1) = 0:1; ~
F (2; 2) = 0:5;
ffi 0 = 0:7; jq 1 (t)j 6 0:1; jq 2 (t)j 6 0:1
N = 1; S =
` 4 0
0 4
'
; G =
` 1 0
0 1
'
; P 0 =
` 1 0
0 1
'
Fig. 2
The worst perturbation p 0 = \Gamma0:7 corresponds to the case when the current coil is located in such a manner
that its magnetic torque is directed along the x 3 ­axis. It follows from Figure 2 that the position when the frame
is perpendicular to the x 3 ­axis provides the best stabilization of pitching oscillations.

Mathematical Modeling of Complex Information Processing Systems 69
The behavior of mean value of the satellite's rotation angle e x1 and its angular velocity e x2 as well as the
worst perturbation q 0
1 (t) for the worst p 0 = \Gamma0:7 are shown in Figure 3.
At the second stage of testing, we perform statistical tests for the stabilization of oscillations with the use
of the control signals u r = \Gamma sign
\Theta (b 00 + b 12 p) —
ff
\Lambda (this control is examined in [8]).
As for the sample mean value of performance criterion (2) and its terminal part corresponding to the
control u r , the dependence of this value on the test number k is represented in Figure 4. As is seen from this
figure, in order to obtain a stable statistical estimate of the functional, we should carry out approximately 60--80
tests. Statistical tests of the control algorithm u 0 = K(p 0 )x obtained at the first stage were also performed.
The time dependence of elements of the covariance matrix of the state vector P 0
ij for the optimal control u 0
(the thick lines) and the covariance values P r
ij for the control u r to be tested (the thin lines) are represented in
Figure 5.
The results of testing are given in Table 1. The sample mean values m J , the sample variances oe J , and
the corresponding confidence intervals I m , I oe (with the confidence probability p ff = 0:95) are shown in the first
row for the best accuracy J 0
0 of criterion (2) (the control u 0 ). The second row contains the results of statistical
tests for the quality functional J 0
r corresponding to the control u r to be tested. The last two rows of the table
contain similar values for the terminal part of functional (2).
Fig. 3 Fig. 4
Fig. 5

70 Mathematical Modeling of Complex Information Processing Systems
Table 1
J m J oe J I m (0.05) I oe (0.05)
J 0
0 19.8 1.037 19.63--20.01 0.92--1.188
J 0
r 29.84 2.05 29.47--30.21 1.82--2.348
J t
0 0357 0.092 0.340--0.374 0.85--1.102
J t
r 3.27 0.96 3.09--3.44 0.85--1.102
Table 2
J m J oe J I m (0.05) I oe (0.05)
J 0
0 40.55 2.22 40.15--40.96 1.97--2.55
J 0
r 66.06 13.05 63.7--68.42 11.58--14.94
J t
0 2.28 0.35 2.21--2.35 0.316--0.407
J t
r 3.01 1.48 2.75--3.28 1.313--1.695
It follows from Table 1 that in terms of the chosen criterion the quality of the control u r is Ÿ = 29:84
19:8 = 1:507
times worse than it is for the ``best'' result of control; with a probability of 0.95, Ÿ belongs to the range
1:472 ! Ÿ ! 1:538.
If we consider only the terminal part of criterion (2), then the results are essentially worse: Ÿ = 9:3. This
is due to the fact that the control u 0 responds to small deviations of the satellite from its stationary position
better than the discrete control u r . If stronger requirements are imposed on the auxiliary control u 0 , then the
results of testing become better. As an example, the results of statistical tests performed with N = 3:5 in
functional (2) are given in Table 2.
As is seen from Table 2, with this performance criterion the result of testing the terminal part of func­
tional (2) is Ÿ = 3:01
2:28 = 1:32 if N = 3:5.
The result of testing for functional (2) (Ÿ = 66:06
40:55 = 1:628) is no better than in the case of the criterion
for N = 1.
REFERENCES
1. V.V. Alexandrov, ``Testing of stabilization quality of nonsteady motions'', Vestn. Mosk. Univ. Matem.
Mekhan., 3: 51--54, 1997.
2. V.V. Alexandrov, S.I. Zlochevskii, S.S. Lemak, and N.A. Parusnikov, Introduction to Dynamics of Control
Systems (in Russian), Moscow, 1993.
3. N.N. Krasovskii, Control of Dynamical Systems (in Russian), Moscow, 1985.
4. M.H.A. Davis, Linear Estimation and Stochastic Control, London, 1977.
5. I.A. Krylov and F.L. Chernous'ko, ``A method of successive approximations for problems of optimalcontrol'',
Zhurn. Vychislit. Matem. i Matem. Fiz., 2, 6: 1132--1138, 1962.
6. V.V. Alexandrov, ``Accumulation of perturbations in linear systems in two coordinates'', Vestn. Mosk.
Univ. Matem. Mekhan., 3: 67--76, 1968.
7. K. Schittkowski, ``The nonlinear programming method of Wilson, Han, and Powell with an augmented
Lagrangian type line search function'', Numer. Math., 38: 115--127, 1981.
8. V.V. Alexandrov and V.N. Zhermolenko, ``Minimax stabilization of a parametrically perturbed oscillatory
system'', Vestn. Mosk. Univ. Matem. Mekhan., 6: 40--43, 1998.