Coordinated Motion Control of Distributed Spacecraft with Relative State Estimation

Wang, Xiaochu; Zhao, Kaichun; You, Zheng

doi:10.1061/(ASCE)AS.1943-5525.0000559

Open access

Technical Papers

Oct 19, 2015

Coordinated Motion Control of Distributed Spacecraft with Relative State Estimation

Authors: Xiaochu Wang, Ph.D. [email protected], Kaichun Zhao, Ph.D. [email protected], and Zheng You, Ph.D. [email protected]Author Affiliations

Publication: Journal of Aerospace Engineering

Volume 29, Issue 3

https://doi.org/10.1061/(ASCE)AS.1943-5525.0000559

PDF

Abstract

This paper addresses the coordination problem of distributed spacecraft systems with time-invariant communication topologies. In order to estimate the neighbor spacecraft’ relative states which are not able to be measured directly, estimating processes are established on each spacecraft by utilizing an unscented Kalman filtering approach with measuring principles of both the pseudorange and the carrier phase. Based on such estimations, a distributed consensus protocol is proposed to address the coordinated motion problem of distributed spacecraft. The effectiveness of the proposed approach is demonstrated through numerical simulations, with a low-Earth-orbit satellites flying scenario.

Introduction

The coordination problem of distributed spacecraft systems has attracted much attention in recent years. Compared to solo systems, additional benefits such as high flexibility and great efficiency can be obtained by having a group of distributed spacecraft work cooperatively.

A fundamental approach to achieve such coordination is consensus, which enables all agents of a system to reach an agreement on certain quantities of interest. Consensus problems have been widely studied in recent years. Jadbabaie et al. (2003) analyzed the coordination behavior of a group of autonomous agents and provided a theoretical explanation for such behavior, and their results indicated that proper control models based only on neighbor information could motivate systems to achieve consensus. Ren et al. (2007) made an overview on information consensus of multiple vehicles by not only reviewing different consensus control laws but also describing several potential specific applications, such as rendezvous, formation stabilization, formation maneuvering, and flocking. Some more challenging results were obtained in Ren (2008) by extending consensus algorithms from first-order models to second-order ones, which made consensus approaches applicable for systems with second-order dynamics. Some practical conditions such as the actuator saturations and the lacks of relative velocity measurements were also discussed. In the context of spacecraft coordination, Ren (2007) considered consensus problems and showed the concept of coordinated motion. By utilizing a consensus control law with local information exchanges, the given formation configuration was maintained during the maneuvers, and all states of spacecraft were assumed available with respect to the Earth Centered Inertial (ECI) frame. In order to obtain accurate formation configurations, orbital reference frames were used to represent relative motions of spacecraft instead of the ECI frame (Chung et al. 2009; Li et al. 2010), and all relative states of spacecraft were also assumed available in this work.

However in practice, different from absolute states (i.e., positions and velocities), relative states (i.e., relative positions and relative velocities) are generally difficult to directly measure due to the limitation of sensors. Therefore, additional work should be done to obtain these relative quantities. Among varieties of methods to do so, although the direct method that simply subtracts between two irrelevant absolute states could be used to calculate the relative state, the results would be quite rough. In contrast, methods based on relative sensing are preferred. Carrier-phase-based differential global positioning system navigation is one such method that is applicable only for missions near the Earth [e.g., Mitchell (2004)]. For more-general missions where the global positioning system may not be available, such as deep-space missions, several alternative sensing technologies have been developed in recent years. Some common such technologies are based on lasers, charge coupled devices, and high-speed cameras. However, laser-based sensing is of high directionality, and image-based sensing requires complex image processing and target identification—with associated failures. Gunnam et al. (2002) developed a new vision-based navigation (VISNAV) sensor at Texas A&M University (College Station, Texas) for relative navigation. VISNAV uses at least four active beacons, which are at known positions in the body frame of one spacecraft, as light sources. It also uses a position-sensitive diode sensor, which is at a known position in the body frame of the other spacecraft, as the light receiver. Thus after gathering all the line-of-sight vectors generated by light sources and receiver, relative states between the two spacecraft can be estimated. Lau et al. (1996) provided an autonomous formation flying sensor (AFF) for relative navigation, which is built around commercially-available transmitter/receiver antennas at known positions in body frames of each spacecraft. By dealing with the measured pseudorange or phase observables, relative states can be solved to extremely high precision by the estimation. In the context of distributed spacecraft systems, since the vision-based method may be blocked by other spacecraft due to the complexity and flexibility of the formations, AFF is considered as the sensor whose sensing range can be from several meters to several kilometers.

This paper considers the coordinated motion problem of distributed spacecraft systems that consist of one free-flying spacecraft and several maneuvering spacecraft. All spacecraft are equipped with the AFF sensors, and besides, only the free-flying spacecraft has the ability to sense its absolute position, while the maneuvering spacecraft do not. The organization of this paper proceeds as follows. First, some background on algebraic graph theory, definition of coordinates, and AFF sensors are given. Then the measurement equations of the AFF sensor, which will be used for filtering, are derived in terms of both pseudorange and phase. Next, the equation of relative motion of spacecraft is reviewed and presented, and a set of unscented kalman filters are designed for relative estimation in a distributed way. During the design, the initial values of filters are given not by a priori assumption but by the solution provided by a nonlinear minimization algorithm. Next, a consensus control law is designed for the orbital coordinative motion. In the simulation, a low-Earth-orbit satellites flying scenario is taken as an instantiation, and subsequent results and analyses are presented. Finally, concluding remarks are made in the last section.

Backgrounds

Algebraic Graph Theory

The bidirectional information exchange among spacecraft can be modeled by undirected graphs. A weighted undirected graph

G = (N, E)

describes the communication links of a finite nonempty set of nodes

N = {1, 2, \dots, n}

with

E \subseteq N \times N

being the edge set. In a weighted undirected graph, an edge

(j, i) \in E

denotes that spacecraft

i

and

j

can obtain information from each other, and they are called neighbors. The weighted adjacency matrix

A = [a_{i j}] \in R^{n \times n}

of a weighted undirected graph is defined such that

a_{i j}

is a positive weight if

(j, i) \in E

, while

a_{i j} = 0

if

(j, i) \notin E

(without loss of generality, it is considered that

a_{i j}

has been turned into the range of

[0, 1]

in this paper). As

(j, i) \in E

implies

(i, j) \in E

, it follows that

a_{i j} = a_{j i}

,

\forall i \neq j

. The Laplacian matrix

L = [l_{i j}] \in R^{n \times n}

associated with the graph is defined such that

l_{i i} = \sum_{j = 1, j \neq i}^{n} a_{i j}

, while

l_{i j} = - a_{i j}

,

\forall i \neq j

. It is clear that

A

is symmetric, and

L

is symmetric, positive, and semidefinite when the undirected graph is connected [see e.g., Ren et al. (2007) for more].

Coordinates

In this paper, the free-flying spacecraft is treated as a reference spacecraft. To describe the absolute position of the reference spacecraft, an ECI frame named J2000 is used, while to describe the relative motion between every two spacecraft, an orbital reference frame named the local vertical/local horizontal (LVLH) frame is used. As shown in Fig. 1,

O_{E} X_{I} Y_{I} Z_{I}

denotes the ECI frame with the origin

O_{E}

attached to the mass center of the Earth;

O_{R} X_{R} Y_{R} Z_{R}

denotes the LVLH frame with the origin

O_{R}

attached to the mass center of the reference spacecraft; Point

O_{i}

denotes the mass center of the

i

th maneuvering spacecraft; and vector

ρ_{i}

denotes the relative position of the

i

th maneuvering spacecraft with respect to the LVLH frame. The LVLH frame is defined such that the

X_{R}

axis is outward along the radial (local vertical),

Y_{R}

is perpendicular to

X_{R}

in the orbit plane in the direction of motion (local horizontal), and

Z_{R}

is along the orbit normal.

AFF Sensor

AFF sensor was first devised for the Deep Space Mission 3 (DS-3) (Lau et al. 1996), where the stabilization of a fixed triangle formation is strictly required. In this mission, since the three spacecraft always face towards each other, AFF is configured with only a

2 π

steradian coverage for the purpose of saving cost and mass. There is also an alternative full configuration for more flexible formations, and in this paper, this full configuration is employed as shown in Fig. 2 for a full

4 π

steradian coverage. For simplicity, the shape of each spacecraft is simplified to a cube, and coordinate

O_{i} X_{i} Y_{i} Z_{i}

is the body fixed frame of spacecraft

i

(

BF - i

). Each spacecraft is equipped with the AFF consisting of eight antennas, of which six are receivers and two are transmitters at known positions in the

BF - i

frame. Symbol

B

is borrowed to indicate the receivers, while Symbol

C

indicates the transmitters. Superscripts are identifiers of distinct transmitters/receivers, and subscripts indicates the owner spacecraft of these antennas (e.g.,

B_{i}^{3}

means Receiver 3 on spacecraft

i

, and

C_{j}^{1}

means the Transmitter 1 on spacecraft

j

). In a single-epoch measurement, ranging codes and carrier phases are generated upon each spacecraft, sent to the neighboring spacecraft, and processed in sequence.

Relative Estimation

Measurement Equations

For convenience, the maneuvering spacecraft is labeled with number

i = 1, 2, \dots, n

in turn, and to achieve a unified format, the reference spacecraft is labeled with number 0. Now consider the transmit–receive relations between two such distinct spacecraft

i

and

j

(

i, j \in {0, 1, 2, \dots, n}

and

i \neq j

) that can communicate with each other. Since there are two transmit antennas mounted on opposing corners on one spacecraft, there will always exist at least one transmit antenna that can be received by the other spacecraft. And in the meanwhile, there will always exist at least three active receive antennas on one spacecraft that can receive signals from the other spacecraft. Therefore, the directed transmit–receive paths between spacecraft

i

and

j

are at least six, among which only six will be selected (can be randomly) to take part in the measurement. Such six transmit–receive paths are illustrated in Fig. 3, where the arrows represent the directions of the information flows.

Fig. 3. Illustration of six transmit–receive paths

Pseudorange measurement equations will now be derived. As shown in Fig. 4, without loss of generality, focusing on one general path where the active receive antenna on spacecraft

i

is denoted as

B_{i}^{r}

and the active transmit antenna on spacecraft

j

is denoted as

C_{j}^{t}

,

r \in {1, 2, \dots, 6}

,

t \in {1, 2}

. According to geometric vector addition, one obtains

{\vec{C_{j}^{t} B_{i}^{r}}}^{I} = {\vec{O_{i} B_{i}^{r}}}^{I} - {\vec{O_{i} O_{j}}}^{I} - {\vec{O_{j} C_{j}^{t}}}^{I}

(1)

where the superscripts of vectors indicate the coordinate frames in which the vectors are expressed, and

I

,

L

,

i

and

j

denote the ECI, LVLH,

BF - i

, and

BF - j

frames, respectively.

Denoting the transformation matrix from frame

F^{1}

to

F^{2}

as

A_{F^{2} F^{1}}

, Eq. (1) can be rewritten as

{\vec{C_{j}^{t} B_{i}^{r}}}^{I} = A_{I i} \cdot {\vec{O_{i} B_{i}^{r}}}^{i} - A_{I L} \cdot {\vec{O_{i} O_{j}}}^{L} - A_{I j} \cdot {\vec{O_{j} C_{j}^{t}}}^{j}

(2)

Then, taking measurement noises and clock errors into consideration, pseudorange measurement equations can be obtained as

‖ {\vec{C_{j}^{t} B_{i}^{r}}}^{I} ‖ = ‖ A_{I i} \cdot {\vec{O_{i} B_{i}^{r}}}^{i} - A_{I L} \cdot {\vec{O_{i} O_{j}}}^{L} - A_{I j} \cdot {\vec{O_{j} C_{j}^{t}}}^{j} ‖ + n_{i j} + ε

(3)

where

‖ \cdot ‖

= 2-norm operator;

n_{i j} = n_{i} - n_{j}

= clock error between spacecraft

i

and

j

, which is converted to the form of distance; and

ε

= measuring noise, which can be modeled to follow the normal distribution. Also note that Eq. (3) represents three equations corresponding to the three paths from spacecraft

j

to

i

.

The paths from spacecraft

i

to

j

are similar, and equations can be obtained by simply switching the variables

i

and

j

in Eq. (3) as

‖ {\vec{C_{i}^{t} B_{j}^{r}}}^{I} ‖ = ‖ A_{I j} \cdot {\vec{O_{j} B_{j}^{r}}}^{j} - A_{I L} \cdot {\vec{O_{i} O_{j}}}^{L} - A_{I i} \cdot {\vec{O_{i} C_{i}^{t}}}^{i} ‖ + n_{j i} + ε

(4)

where

n_{j i} = - n_{i j}

. In Eqs. (3) and (4),

{\vec{O_{i} B_{i}^{r}}}^{i}

,

{\vec{O_{j} C_{j}^{t}}}^{j}

,

{\vec{O_{j} B_{j}^{r}}}^{j}

, and

{\vec{O_{i} C_{i}^{t}}}^{i}

are known vectors,

A_{I i} = A_{i I}^{- 1}

and

A_{I i} = A_{i I}^{- 1}

with

A_{i I}

and

A_{j I}

are known as the attitude matrices of spacecraft

i

and

j

,

A_{I L}

is the known transformation matrix from LVLH frame to ECI frame, while

{\vec{O_{i} O_{j}}}^{L}

and

n_{i j}

are unknown quantities to be solved.

As shown in Fig. 5, the differential phase measurement is abstracted. Denote

B_{i}^{r_{1}}

and

B_{i}^{r_{2}}

as the two receivers on spacecraft

i

that receive signals from the same transmitter

C_{j}^{t}

on spacecraft

j

. Then difference the two received signals, the differential phase

Δ ϕ

can be measured. Define line-of-sight vectors as

l o s_{1} = \frac{\vec{C_{j} B_{i}^{r_{1}}}}{‖ \vec{C_{j} B_{i}^{r_{1}}} ‖}, l o s_{2} = \frac{\vec{C_{j} B_{i}^{r_{2}}}}{‖ \vec{C_{j} B_{i}^{r_{2}}} ‖}

(5)

Then the differential phase measurement equations can be expressed as

Δ ϕ = \vec{B_{i}^{r_{1}} B_{i}^{r_{2}}} \cdot l o s_{2} + β + η - δ

(6)

where

β

= integer ambiguity;

η

= measurement noise; and

δ

= small distance, which can be calculated as

δ = ‖ \vec{C_{j} B_{i}^{r_{1}}} ‖ - ‖ \vec{C_{j} B_{i}^{r_{1}}} ‖ \cdot \cos α = ‖ \vec{C_{j} B_{i}^{r_{1}}} ‖ (1 - l o s_{1} \cdot l o s_{2})

(7)

Since three receive antennas can only make two combinations, there are four differential phase measurements in a single epoch, and hence the equations are four. Combining six pseudorange measurement equations and four differential phase measurement equations yields the final measurement equations.

Relative Motion Equations

As defined previously, the LVLH frame attached to the free flying spacecraft is used as the orbital reference frame in this paper. Although the relative motion can be also described in similar orbital frames attached to any one of the maneuvering spacecraft, the maneuvering spacecraft will be maneuvering in general so each of them would not likely move in a fixed Keplerian orbit. Therefore, it is more accurate to use the LVLH frame attached to the free-flying spacecraft as opposed to the maneuvering spacecraft. Besides, although the relative motion is being studied, the rough orbital absolute information is needed due to the orbital characteristics, so in this paper it is assumed that the free-flying spacecraft is equipped with an absolute position sensor (e.g., the global positioning system) to provide a reference absolute position for the spacecraft cluster.

For simplicity, in this paper, the free-flying spacecraft is assumed to move in a nearly circular orbit. Using Kepler’s third law, the orbit angular velocity can be given by

ω_{o} = \sqrt{\frac{μ_{e}}{R_{o}^{3}}}

(8)

where

ω_{o}

= orbit angular velocity;

μ_{e}

= geocentric gravitational constant (398,

600.4418 {km}^{3} s^{- 2}

); and

R_{o}

= orbit radius. Denote the relative position of the

i

th maneuvering spacecraft with respect to the reference spacecraft in the LVLH frame as

ρ_{i} = {[x_{i}, y_{i}, z_{i}]}^{T}

. When the relative position between a maneuvering spacecraft

i

and the reference spacecraft are small compared to the orbit radius, the relative dynamics, known as the CW equations, can be derived as [see, Clohessy and Wiltshire (1960) and Lee et al. (2014) for more]

{\ddot{x}}_{i} - 2 ω_{o} {\dot{y}}_{i} - 3 ω_{o}^{2} x_{i} = \frac{u_{i x} + d_{i x}}{m_{i}} {\ddot{y}}_{i} + 2 ω_{o} {\dot{x}}_{i} = \frac{u_{i y} + d_{i y}}{m_{i}} {\ddot{z}}_{i} + ω_{o}^{2} z_{i} = \frac{u_{i z} + d_{i z}}{m_{i}}

(9)

where

m_{i}

= mass of the

i

th maneuvering spacecraft;

u_{i x}

,

u_{i y}

, and

u_{i z}

= control forces; and

d_{i x}

,

d_{i y}

, and

d_{i z}

= disturbance forces. Although Eq. (9) is an approximation and is not applicable for long-term orbital motion propagation, it can be used to predict the relative state in the next moment with enough precision. Furthermore, motion within the orbital plane (i.e.,

x

and

y

axes) are coupled, while the motion perpendicular to the orbital plane (i.e.,

z

axis) is independent from that within the orbital plane, and actually it is periodic vibratory. Denote

v_{i} = {\dot{ρ}}_{i} = {[{\dot{x}}_{i}, {\dot{y}}_{i}, {\dot{z}}_{i}]}^{T}

and define

X_{i} = {[ρ_{i}^{T}, v_{i}^{T}, n_{i}^{T}]}^{T}

, then the equation of

X_{i}

can be rewritten in the matrix form as

{\dot{X}}_{i} = F X_{i} + G U_{i} + G D_{i}

(10)

with Fig. 6

G = [\begin{matrix} 0_{3 \times 3} \\ I_{3 \times 3} \\ 0_{1 \times 3} \end{matrix}]

(11)

U_{i} = \frac{1}{m_{i}} {[u_{i x}, u_{i y}, u_{i z}]}^{T}, D_{i} = \frac{1}{m_{i}} {[d_{i x}, d_{i y}, d_{i z}]}^{T}

(12)

where

0_{3 \times 3}

is a

3 \times 3

zero matrix, and

I_{3 \times 3}

is a

3 \times 3

identity matrix. Furthermore, by defining

X_{i j} = X_{i} - X_{j}

,

U_{i j} = U_{i} - U_{j}

, and

D_{i j} = D_{i} - D_{j}

, the relative equations between two distinct maneuvering spacecraft

i

and

j

(

i, j \in {1, 2, \dots, n}

and

i \neq j

) can be written as

{\dot{X}}_{i j} = F X_{i j} + G U_{i j} + G D_{i j}

(13)

Fig. 6. The system matrix F in the state differential equation Eq. (10)

Estimation

In this subsection, the method for relative estimations among multiple spacecraft will be presented. First, how the estimation process works at the system level will be shown. Since the spacecraft system is distributed and a host spacecraft does not exist, the estimation could be processed by many small estimators working together in a distributed way. It is clear that the number of small estimators is equal to that of undirected edges in the topology. To describe the allocation of these estimators, the sign matrix

S = [s_{i j}] \in R^{n \times n}

is defined such that

s_{i j} = 1

if an estimator for

X_{i j}

is running on spacecraft

i

, while

s_{i j} = 0

otherwise. Define

s_{i} = \sum_{k = 0}^{n} s_{i k}

, then the problem of allocating the estimators evenly is equivalent to the problem of trying to average

s_{i}

,

i = 0, 1, 2, \dots, n

. Hence a method to get an averaging decision is proposed as follows:

1.

For each spacecraft

i \in {0, 1, 2, \dots, n}

, initialize

s_{i j} = 0

,

j = 0, 1, 2, \dots, n

;

2.

During one communication, for each pair of neighbor spacecraft

i

and

j

, a random decision is made such that either

s_{i j} = 1

,

s_{j i} = 0

or

s_{i j} = 0

,

s_{j i} = 1

;

3.

During a next communication, for each pair of neighbor spacecraft

i

and

j

, an improved decision is made by

s_{i j} (k) = s_{i j} (k - 1) - \frac{1}{n + 1} a_{i j} [s_{i} (k - 1) - s_{j} (k - 1)] s_{i j} (k) = 1, if s_{i j} (k) > 1 s_{i j} (k) = 0, if s_{i j} (k) < 0

(14)

where (

k

) = moment;

4.

Repeat Step. 3 for a given while. And then by rounding

s_{i j}

to the nearest integers, one obtains the final decision.

It can be seen from Eq. (14) that if spacecrafts

i

and

j

are neighbors, there comes

a_{i j} > 0

, so that when

s_{i} > s_{j}

(which means the work of spacecraft

i

is heavier than that of spacecraft

j

),

s_{i j}

will decrease (meanwhile

s_{j i}

will increase) to offload a little work from

i

to

j

, and when

s_{i} < s_{j}

, then

s_{i j}

will increase to help share some work from

j

to

i

. Thus after some time, the amount of computations would be well offloaded to every spacecraft. Eq. (14) also constrains that

0 \leq s_{i j} \leq 1

, so after rounding towards the nearest integers, the final decision can be obtained. The gain “

1 / n + 1

” was chosen to guarantee the convergence of the averaging iteration. Actually, by denoting the gain as

κ

, then the gain can be chosen as any values satisfying

0 < κ \leq \min_{0 \leq i \leq n} (1 / \sum_{j = 0}^{n} a_{i j})

. This choosing principle can be proved as follows. Replace the gain with

κ

, and the first equation of Eq. (14) can be written as

s_{i j} (k) = s_{i j} (k - 1) - a_{i j} \cdot κ \cdot [s_{i} (k - 1) - s_{j} (k - 1)]

(15)

One can obtain

\sum_{j = 0}^{n} s_{i j} (k) = \sum_{j = 0}^{n} {s_{i j} (k - 1) - a_{i j} \cdot κ \cdot [s_{i} (k - 1) - s_{j} (k - 1)]}

by summing Eq. (15) together. Integrate these equations and get the matrix form as

s_{k} = (I - κ \cdot L) \cdot s_{k - 1}

(16)

where

s_{k} = {[s_{0} (k), s_{1} (k), \dots, s_{n} (k)]}^{T}

,

I

is the

n + 1

-order identity matrix. Denote by

λ_{i}

the

i

th eigenvalue of

(κ \cdot L)

in an increasing sequence and by

μ_{i}

the

i

th eigenvalue of

Θ ≜ (I - κ \cdot L)

, and it is easy to know that

μ_{i} = 1 - λ_{i}

. Choose the gain following the principle of

0 < κ \leq 1 / \sum_{j = 0}^{n} a_{i j}

and the graph is connected; then it follows from Ren et al. (2007) that there is a single 0 in the eigenvalues of

(κ \cdot L)

while all the other eigenvalues are positive but no more than 1, or in formula form,

0 = λ_{0} < λ_{1} \leq \dots \leq λ_{n} \leq 1

. By doing so, the eigenvalues of

Θ

satisfy

0 \leq μ_{n} \leq \dots \leq μ_{1} < μ_{0} = 1

, which implies the spectral radius of

Θ

is 1. It follows from Horn and Johnson (2012) that,

\lim_{k \to \infty} Θ^{k} = 1 \cdot ξ^{T}

, where

ξ^{T}

is the left eigenvector of

Θ

corresponding to the eigenvalue

μ_{0}

. Therefore after enough loops,

s_{k}

converges following

\lim_{k \to \infty} s_{k} = \lim_{k \to \infty} Θ^{k} s_{0} = 1 \cdot ξ^{T} \cdot s_{0} = 1 \cdot \bar{s}

, where

\bar{s}

denotes a constant. It means

s_{1} = s_{2} = \dots = s_{n}

after convergence, which indicates that all the spacecraft finally reach consensus on the division of work.

After offloading computations to every spacecraft by the averaging decision, each spacecraft will then activate a corresponding estimator according to

s_{i j}

. These estimators are designed based on unscented Kalman filters as follows:

1.

Initialization: At moment 0, initialize the estimation of relative state

{\hat{X}}_{i j} (0)

, process noise matrix

Q_{i j} (0)

, measurement noise matrix

R_{i j} (0)

, and error covariance matrix

P_{i j} (0)

. Further analyses and principles that accomplish the initialization will be given following the design.

2.

Sampling weights assignment: The weight assignment follows the minimal skew sampling method (Julier and Uhlmann 2002). Since

X_{i j}

contains seven quantities, the number of sigma points should be

7 + 2 = 9

. The weights are chosen as

W_{m} = {\begin{cases} \frac{1 - W_{0}}{27}, & m = 1, 2 \\ 2^{m - 1} W_{1}, & m = 3, \dots, 8 \end{cases}

(17)

with

W_{0}

being a initial nonnegative weight, and in this paper

W_{0} = 0.6

.

3.

Unbiased sigma points construction of u: Since the covariance is already reflected in matrices

P_{i j}

and

Q_{i j}

, the sigma points should be unbiased. Denote the nine sigma points as

S_{m}^{7}

,

m = 0, \dots, 8

respectively, then their values can be given by following iterations

S_{i}^{j} = {\begin{cases} [\begin{matrix} S_{0}^{j - 1} \\ 0 \end{matrix}], & i = 0 \\ [\begin{matrix} S_{i}^{j - 1} \\ - \frac{1}{\sqrt{2 W_{j + 1}}} \end{matrix}], & i = 1, \dots, j \\ [\begin{matrix} 0_{(j - 1) \times 1} \\ \frac{1}{\sqrt{2 W_{j + 1}}} \end{matrix}], & i = j + 1 \end{cases}

(18)

with three priori values

S_{01} = 0, S_{11} = - \frac{1}{\sqrt{2 W_{1}}}, S_{21} = \frac{1}{\sqrt{2 W_{1}}}

(19)

4.

Global propagation: At moment

k - 1

, the purpose of the global propagation is to predict

{\hat{X}}_{i j} (k)

roughly based on

{\hat{X}}_{i j} (k - 1)

, and the prediction is denoted as

{\hat{X}}_{i j} (k, k - 1)

. Then following Eq. (13), the prediction is given by

{\hat{X}}_{i j} (k, k - 1) = Φ \cdot {\hat{X}}_{i j} (k - 1) + Δ T \cdot {G U}_{i j}

(20)

where

Δ T

= time interval between estimations; and

Φ

= state transition matrix. In practice, it is of sufficient precision to preserve the first moment only, and hence one obtains

Φ \approx I_{7 \times 7} + F \cdot Δ T

(21)

5.

Error propagation: The purpose of error propagation is to evaluate the errors caused by uncertainties of

{\hat{X}}_{i j} (k - 1)

during the global propagation. Since Eq. (20) is linear, the prediction itself can be omitted in this process. Based on the unbiased sigma points, perturbations before global propagation can be constructed as

χ_{m} (k - 1) = \sqrt{P_{i j} (k - 1) + Q_{i j} (k - 1)} \cdot S_{m}^{7}, m = 0, \dots, 8

(22)

which will generate propagation errors as

χ_{m} (k, k - 1) = Φ \cdot χ_{m} (k - 1), m = 0, \dots, 8

(23)

6.

Measurement prediction: The measurement equations Eqs. (3), (4), and (6) are all functions of relative state

\hat{X}

, hence for clarity, they are treated as a whole, which is referred to as the following function

y = f (\hat{X})

(24)

Therefore, the predictions of measurements are given by

Z_{m} (k, k - 1) = f [{\hat{X}}_{i j} (k, k - 1) + χ_{m} (k, k - 1)]

(25)

7.

Covariance prediction: Covariances are predicted as

P_{i j} (k, k - 1) = \sum_{m = 0}^{8} W_{m} [χ_{m} (k, k - 1) - δ {\hat{X}}_{i j} (k, k - 1)] \cdot {[χ_{m} (k, k - 1) - δ {\hat{X}}_{i j} (k, k - 1)]}^{T} + Q_{i j} (k - 1) P_{i j}^{y y} (k, k - 1) = \sum_{m = 0}^{8} W_{m} [Z_{m} (k, k - 1) - Z_{i j} (k, k - 1)] \cdot {[Z_{m} (k, k - 1) - Z_{i j} (k, k - 1)]}^{T} P_{i j}^{x y} (k, k - 1) = \sum_{m = 0}^{8} W_{m} [χ_{m} (k, k - 1) - δ {\hat{X}}_{i j} (k, k - 1)] \cdot {[Z_{m} (k, k - 1) - Z_{i j} (k, k - 1)]}^{T}

(26)

where

δ {\hat{X}}_{i j} (k, k - 1) = \sum_{m = 0}^{8} W_{m} χ_{m} (k, k - 1) Z_{i j} (k, k - 1) = \sum_{m = 0}^{8} W_{m} Z_{m} (k, k - 1)

(27)

8.

Update: Once the true measurements

Z_{i j} (k)

occur at moment

k

, the estimation will be updated as

{\hat{X}}_{i j} (k) = {\hat{X}}_{i j} (k, k - 1) + K_{i j} [Z_{i j} (k) - Z_{i j} (k, k - 1)]

(28)

where

K_{i j} = P_{i j}^{x y} (k, k - 1) {[P_{i j}^{y y} (k, k - 1) + R_{i j} (k)]}^{- 1}

(29)

9.

Assignment for next moment: After the estimation at moment

k

, the following assignments will be done to get ready for the estimation at moment

k + 1

, which are

P_{i j} (k) = P_{i j} (k, k - 1) - K_{i j} [P_{i j}^{y y} (k, k - 1) + R_{i j} (k)] K_{i j}^{T} Q_{i j} (k) = Q_{i j} (k - 1) R_{i j} (k) = R_{i j} (k - 1)

(30)

The initializations of estimators are assumed to be priori known in general studies, however, this assumption is difficult to realize in practice. In this paper, practical principles of initialization will be presented.

{\hat{X}}_{i j}

contains three quantities which are

{\hat{ρ}}_{i j}

,

{\hat{v}}_{i j}

, and

{\hat{n}}_{i j}

respectively. Regarding

{\hat{n}}_{i j}

, the work done by Lau et al. (1996) show that clock error could be measured by a built-in approach with the precision of 0.1 ns, which is equivalent 0.03 m after being multiplied by the speed of light. Therefore,

{\hat{n}}_{i j} = 0

can be initialized with the bound

σ_{n} = 0.03 m

. Regarding

{\hat{v}}_{i j}

, the cluster of distributed spacecraft are assumed have been already launched approximately to a common orbit such that they can communicate with neighbors, and hence, their initial relative velocities are bounded. The bound can be approximately estimated simply. Since

‖ v_{i j} ‖ \approx (‖ ρ_{i j} ‖ / R_{o}) V_{o}

with

R_{o}

bigger than the Earth’s radius,

V_{o}

smaller than the first cosmic velocity, and

‖ ρ_{i j} ‖

no bigger than several kilometers (in fact, the maximum effective distance of AFF is commonly 2 km) in general, there comes the inequality that

‖ v_{i j} ‖ < (2 / 6400) \times 7900 < 3 m / s

. Therefore,

{\hat{v}}_{i j} = {[0, 0, 0]}^{T}

can be initialized with the bound

σ_{v} = 3 m / s

. Regarding

{\hat{ρ}}_{i j}

, it can be solved roughly based on Eqs. (3) and (4) ignoring the small

n_{i j}

. Thus the initialization of

{\hat{ρ}}_{i j}

becomes the optimization problem optimizing

ρ_{i j}

to make the results of Eqs. (3) and (4) best fit the pseudorange measurements. Since the number of equations (i.e., 6) is more than the number of variables (i.e., 3), this problem is a determined system. Therefore, many nonlinear minimization methods can be used for help. In this paper, the trust-region method is used to make a large-scale optimization, and then

{\hat{ρ}}_{i j}

can be initialized to the optimization result with a confidence bound

σ_{ρ}

. Detailed explanations about the algorithm procedure are omitted here due to its universal applicability, which can be found in many good articles, e.g., Coleman and Li (1996) and Conn et al. (2000). So

{\hat{X}}_{i j} (0)

can be initialized until now, and subsequently

P_{i j} (0)

,

Q_{i j} (0)

and

R_{i j} (0)

will be initialized as

P_{i j} (0) = 2 [\begin{matrix} σ_{ρ}^{2} I_{3 \times 3} \\ σ_{v}^{2} I_{3 \times 3} \\ σ_{n}^{2} \end{matrix}] Q_{i j} (0) = 2 Δ T \cdot G \cdot [\begin{matrix} σ_{d x}^{2} \\ σ_{d y}^{2} \\ σ_{d z}^{2} \end{matrix}] \cdot G^{T} R_{i j} (0) = 2 [\begin{matrix} σ_{ε}^{2} I_{6 \times 6} \\ σ_{η}^{2} I_{4 \times 4} \end{matrix}]

(31)

where

σ_{d x}

,

σ_{d y}

, and

σ_{d z}

= confidence bound of the disturbance forces respectively;

σ_{ε}

= confidence bound of pseudorange measurement noises;

σ_{η}

= confidence bound of differential phase measurement noises; and the constant factor 2 here means that corresponding covariances are doubled when differencing

X_{i}

and

X_{j}

.

Distributed Control Based on Estimations

Analyses on Design of Control Law

In this section, a distributed control law will be proposed for relative orbital motion. Suppose a formation configuration is desired for the spacecraft cluster, and it is described by a set of

θ_{i j}

,

i

,

j \in {0, 1, 2, \dots, n}

, which is to say,

θ_{i j}

denotes the desired relative position between spacecraft

i

and

j

. Disturbance forces are generally considered as external disturbances, and hence the system dynamics used in design of the control law can be obtained from the governing dynamic Eq. (9) by ignoring the disturbance terms as

{\dot{ρ}}_{i} = v_{i} {\dot{v}}_{i} = U_{i} + M_{1} ρ_{i} + M_{2} v_{i}, i = 1, \dots, n

(32)

where

M_{1} = [\begin{matrix} 3 ω_{o}^{2} & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & - ω_{o}^{2} \end{matrix}], M_{2} = [\begin{matrix} 0 & 2 ω_{o} & 0 \\ - 2 ω_{o} & 0 & 0 \\ 0 & 0 & 0 \end{matrix}]

(33)

It can be seen that the dynamics presented by Eq. (32) contain the two orbital terms

M_{1} ρ_{i}

and

M_{2} v_{i}

, which are generated by uncontrolled orbital motion. And hence, he effects caused should be eliminated by these terms as much as possible. However, although the differential relative states

ρ_{i j}

and

v_{i j}

can be obtained by estimation,

ρ_{i}

and

v_{i}

themselves are not direct available to spacecraft

i

due to the fact spacecraft

i

will only communicate with its neighbors and could not know the global configuration. Hence, prior to the design of control law, a distributed method is proposed, which is inspired by Cao et al. (2010) and Meng et al. (2010), to help in obtaining the estimation of

ρ_{i}

.

Observing $ρ_{i}$ Based on $ρ_{i j}$

The observing method is proposed as

{\dot{\hat{r}}}_{i} = - γ_{r} sgn [\sum_{j = 1}^{n} a_{i j} ({\hat{r}}_{i} - {\hat{r}}_{j} - {\hat{ρ}}_{i j})], i = 1, \dots, n

(34)

where

sgn (\cdot)

= signum function;

{\hat{r}}_{i}

= estimation of

ρ_{i}

;

{\hat{ρ}}_{i j}

= estimation of

ρ_{i j}

and is already known after the process shown in the previous section; and

γ_{r}

= positive gain satisfying

γ_{r} > \max_{1 \leq i \leq n} {‖ v_{i} ‖}_{\infty}

. The theoretical proof of the method is as follows. Denote

{\tilde{r}}_{i} = {\hat{r}}_{i} - ρ_{i}

,

i = 1, \dots, n

, and then Eq. (34) can be rewritten as

{\dot{\tilde{r}}}_{i} = - {\dot{ρ}}_{i} - γ_{r} sgn [\sum_{j = 1}^{n} a_{i j} ({\tilde{r}}_{i} - {\tilde{r}}_{j} - {\hat{ρ}}_{i j} + ρ_{i} - ρ_{j})] = - {\dot{ρ}}_{i} - γ_{r} sgn [\sum_{j = 1}^{n} a_{i j} ({\tilde{r}}_{i} - {\tilde{r}}_{j})], i = 1, \dots, n

(35)

where the fact that

{\hat{ρ}}_{i j} = ρ_{i j} = ρ_{i} - ρ_{j}

is used to develop the equality. Defining

\tilde{r} = {[{\tilde{r}}_{1}^{T}, {\tilde{r}}_{2}^{T}, \dots, {\tilde{r}}_{n}^{T}]}^{T}

and

ρ = {[ρ_{1}^{T}, ρ_{2}^{T}, \dots, ρ_{n}^{T}]}^{T}

, then Eq. (35) can be expressed in the matrix form as

\dot{\tilde{r}} = - γ_{r} sgn [(L \otimes I_{3 \times 3}) \tilde{r}] - \dot{ρ}

(36)

where

\otimes

is the Kronecker product operator. The subsequent proof can follow the process given by Meng et al. (2010). Since

L

is positive and semidefinite, consider the Lyapunov function candidate

V = \frac{1}{2} {\tilde{r}}^{T} (L \otimes I_{3 \times 3}) \tilde{r}

(37)

Then the derivation of

V

along the trajectories of Eq. (36) is

\dot{V} = {\tilde{r}}^{T} (L \otimes I_{3 \times 3}) \dot{\tilde{r}} = {\tilde{r}}^{T} (L \otimes I_{3 \times 3}) {- γ_{r} sgn [(L \otimes I_{3 \times 3}) \tilde{r}] - \dot{ρ}} \leq - γ_{r} {‖ (L \otimes I_{3 \times 3}) \tilde{r} ‖}_{1} + {‖ \dot{ρ} ‖}_{\infty} \cdot {‖ (L \otimes I_{3 \times 3}) \tilde{r} ‖}_{1} = - (γ_{r} - {‖ \dot{ρ} ‖}_{\infty}) \cdot {‖ (L \otimes I_{3 \times 3}) \tilde{r} ‖}_{1} \leq 0

(38)

where the Holders inequality is used to obtain the first Holder’s inequality and the fact

γ_{r} > \max_{1 \leq i \leq n} {‖ v_{i} ‖}_{\infty} = {‖ \dot{ρ} ‖}_{\infty}

is used to obtain the last inequality. It then follows from Lyapunov stability theory that Eq. (35) is asymptotically stable and

\tilde{r} \to 0

as

t \to \infty

, and by Lasalle’s Invariance principle, it follows that

{\hat{r}}_{i} \to ρ_{i}

,

i = 1, \dots, n

as

t \to \infty

. So there must exist such a time

\bar{t}

that satisfies

{\hat{r}}_{i} = ρ_{i}

,

i = 1, \dots, n

with the margin of errors very small for all

t \geq \bar{t}

. In fact, it can be concluded from Cao and Ren (2012) that

\bar{t} < \sqrt{{\tilde{r}}^{T} (0) (L \otimes I_{3 \times 3}) \tilde{r} (0)} \sqrt{λ_{\max} (L)} / (γ_{r} - {‖ \dot{ρ} ‖}_{\infty}) λ_{\min} (L)

.

Design of Control Law

The control law is proposed as follows:

U_{i} = {\begin{cases} 0_{3 \times 1}, & t \leq \bar{t} \\ - \sum_{j = 0}^{n} a_{i j} γ_{ρ} \tanh ({\hat{ρ}}_{i j} - θ_{i j}) - \sum_{j = 0}^{n} a_{i j} γ_{v} \tanh ({\hat{v}}_{i j}) - M_{1} {\hat{r}}_{i}, & t > \bar{t} \end{cases}, i = 1, \dots, n

(39)

where

γ_{ρ}

= positive gain of the position term;

γ_{v}

= positive gain of the velocity term; and

\tanh (\cdot)

= hyperbolic tangent function and is used here to limit the output since the control forces are bounded in practice [e.g., Meng et al. (2013) and Yang et al. (2014)]. Using control law Eq. (40) for system dynamics Eq. (33), then the system is asymptotically stable. The proof is shown in the following. For

t \leq \bar{t}

, it is clear that system is always bounded, and for

t > \bar{t}

, by combining Eqs. (32) and (39), one obtains

{\dot{ρ}}_{i} = v_{i} {\dot{v}}_{i} = - \sum_{j = 0}^{n} a_{i j} γ_{ρ} \tanh ({\hat{ρ}}_{i j} - θ_{i j}) - \sum_{j = 0}^{n} a_{i j} γ_{v} \tanh ({\hat{v}}_{i j}) + M_{2} v_{i}, i = 1, \dots, n

(40)

where term

M_{1} {\hat{r}}_{i}

counteracts the term

M_{1} ρ_{i}

since

{\hat{r}}_{i} = ρ_{i}

for

t \geq \bar{t}

. Consider the Lyapunov function candidate

V = \frac{1}{2} \sum_{i = 1}^{n} \sum_{j = 1}^{n} {{a_{i j} γ_{ρ} 1_{3}}^{T} \ln [\cosh (ρ_{i j} - θ_{i j})]} + 12 \sum_{i = 1}^{n} v_{i}^{T} v_{i} + \sum_{i = 1}^{n} {{a_{i 0} γ_{ρ} 1_{3}}^{T} \ln [\cosh (ρ_{i 0} - θ_{i 0})]}

(41)

where

\ln (\cdot)

= natural logarithm function; and

\cosh (\cdot)

= hyperbolic cosine function. Then the derivation of

V

along the trajectories of Eq. (40) is

\dot{V} = \frac{1}{2} \sum_{i = 1}^{n} \sum_{j = 1}^{n} [a_{i j} γ_{ρ} v_{i j}^{T} \tanh (ρ_{i j} - θ_{i j})] + \sum_{i = 1}^{n} v_{i}^{T} [- \sum_{j = 0}^{n} a_{i j} γ_{ρ} \tanh (ρ_{i j} - θ_{i j}) - \sum_{j = 0}^{n} a_{i j} γ_{v} \tanh (v_{i j}) + M_{2} v_{i}] + \sum_{i = 1}^{n} [a_{i 0} γ_{ρ} v_{i 0}^{T} \tanh (ρ_{i 0} - θ_{i 0})] = \frac{1}{2} \sum_{i = 1}^{n} \sum_{j = 1}^{n} [a_{i j} γ_{ρ} v_{i j}^{T} \tanh (ρ_{i j} - θ_{i j})] - \sum_{i = 1}^{n} \sum_{j = 1}^{n} [a_{i j} γ_{ρ} v_{i}^{T} \tanh (ρ_{i j} - θ_{i j})] - \sum_{i = 1}^{n} [v_{i}^{T} a_{i 0} γ_{ρ} \tanh (ρ_{i 0} - θ_{i 0})] - \sum_{i = 1}^{n} \sum_{j = 1}^{n} [a_{i j} γ_{ρ} v_{i}^{T} \tanh (v_{i j})] - \sum_{i = 1}^{n} [v_{i}^{T} a_{i 0} \tanh (v_{i})] + \sum_{i = 1}^{n} (v_{i}^{T} M_{2} v_{i}) + \sum_{i = 1}^{n} [a_{i 0} γ_{ρ} v_{i 0}^{T} \tanh (ρ_{i 0} - θ_{i 0})]

(42)

By switching the dummy variables

i

and

j

in the summation function, one obtains

\frac{1}{2} \sum_{i = 1}^{n} \sum_{j = 1}^{n} [a_{i j} γ_{ρ} v_{i j}^{T} \tanh (ρ_{i j} - θ_{i j})] = \sum_{i = 1}^{n} \sum_{j = 1}^{n} [a_{i j} γ_{ρ} v_{i}^{T} \tanh (ρ_{i j} - θ_{i j})]

(43)

\sum_{i = 1}^{n} \sum_{j = 1}^{n} [a_{i j} γ_{ρ} v_{i}^{T} \tanh (v_{i j})] = \frac{1}{2} \sum_{i = 1}^{n} \sum_{j = 1}^{n} [a_{i j} γ_{ρ} v_{i j}^{T} \tanh (v_{i j})]

(44)

also, it can be seen from the definition that

v_{i 0} \equiv v_{i}

, then Eq. (42) can be developed as

\dot{V} = - \frac{1}{2} \sum_{i = 1}^{n} \sum_{j = 1}^{n} [a_{i j} γ_{ρ} v_{i j}^{T} \tanh (v_{i j})] - \sum_{i = 1}^{n} [v_{i}^{T} a_{i 0} \tanh (v_{i})] + \sum_{i = 1}^{n} (v_{i}^{T} M_{2} v_{i})

(45)

As

v_{i}^{T} M_{2} v_{i} \equiv 0

and

\tanh (\cdot)

are odd functions, and then one can obtain

\dot{V} = - \frac{1}{2} \sum_{i = 1}^{n} \sum_{j = 1}^{n} [a_{i j} γ_{ρ} v_{i j}^{T} \tanh (v_{i j})] - \sum_{i = 1}^{n} [v_{i}^{T} a_{i 0} \tanh (v_{i})] \leq 0

(46)

So the system is asymptotically stable and the equilibrium is achieved when

ρ_{i j} = θ_{i j}

and

v_{i j} = 0_{3 \times 1}

,

i = 1, \dots, n

as

t \to \infty

, and thus the formation configuration of the spacecraft system will be formed.

Simulation

A spacecraft system that contains one free-flying spacecraft and four maneuvering spacecraft is taken as an example. A general communication topology is set for the system as shown in Fig. 7, which is fixed, undirected, and connected with communication link edges weighted by 0.8. The five spacecraft are set such that they move in a line initially in a same nearly-circular orbit with orbital semi-major axis 7378.137 km, eccentricity

2.41881 \times 10^{- 16}

, inclination 60°, and right ascension of ascending node 0°. The desired formation offsets and other parameters are given as shown in Table 1. Denote the time when estimation is started as time origin, then the simulation is scheduled as follows. The establishment of communication links and allocation of estimators are processed during

t \in [- 60, 0] s

, all estimators and all observing processes for rough

ρ_{i}

work during

t \in [0, 6000] s

, and the distributed control law works during

t \in [60, 6000] s

. Simulation results are as follows.

Fig. 7. Communication topology of the system

Table 1. Simulation Parameters

Parameter	Value
Desired offset between Spacecrafts 0 and 1, $θ_{10}$	$θ_{10} = {[300, 1000, 0]}^{T} m$
Desired offset between Spacecrafts 1 and 2, $θ_{21}$	$θ_{21} = {[200, 0, 0]}^{T} m$
Desired offset between Spacecrafts 2 and 3, $θ_{32}$	$θ_{32} = {[0, - 200, 0]}^{T} m$
Desired offset between Spacecrafts 3 and 4, $θ_{43}$	$θ_{43} = {[- 200, 0, 0]}^{T} m$
Communication and estimation time interval, $Δ T$	$Δ T = 1 s$
Spacecraft masses, $m_{i}$	$m_{i} = 100 kg$ , $i = 0, \dots, 5$
Disturbance forces, $σ_{x}$ , $σ_{y}$ , $σ_{z}$	$σ_{x} = σ_{y} = σ_{z} {= 10}^{- 6} m / s^{2}$
Pseudorange measurement noises, $σ_{ε}$	$3.33 \times 10^{- 3} m$
Differential phase measurement noises, $σ_{η}$	$σ_{η} = 3.33 \times 10^{- 6} m$
Initial estimation errors of position, $σ_{ρ}$	$σ_{ρ} = 40 m$
Initial estimation errors of position, $σ_{v}$	$σ_{v} = 3 m / s$
Initial estimation clock errors, $σ_{n}$	$σ_{n} = 3 \times 10^{- 2} m$
Gains in Eq. (34), $γ_{r}$	$γ_{r} = \max {1, (600 e^{- 0.03 t} - 200)}$
Gains in Eq. (39), $γ_{ρ}$ , $γ_{v}$	$γ_{ρ} = 0.06$ , $γ_{v} = 0.3$

At time

t = - 60 s

, a random decision is made such that

S = [\begin{matrix} 0 & 0 & 0 & 0 & 0 \\ 1 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 1 \\ 0 & 0 & 0 & 0 & 0 \end{matrix}]

. Following the definition of matrix

S

, it can be seen that Spacecrafts 1, 3, and 5 are idle, while Spacecrafts 2 and 4 both carry out two estimators. At time

t = 0 s

after the allocation, final decision is obtained and

S = [\begin{matrix} 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 \end{matrix}]

, which shows that the four estimators are successfully distributed evenly with Spacecraft 0 estimating

X_{01}

, Spacecraft 1 estimating

X_{12}

, Spacecraft 3 estimating

X_{32}

, and Spacecraft 4 estimating

X_{43}

.

Fig. 8 shows the 2D projections of trajectories of the four maneuvering spacecraft onto the

x - y

plane of the LVLH frame, where the points of the trajectories represent the mass centers of the spacecraft. Fig. 9 shows the detailed 3D position data. From both Figs. 8 and 9, it can be seen that the four maneuvering spacecraft have been successfully controlled to the desired formation configuration. Furthermore, the configuration errors are no more than 1 cm. Controlling is based on estimation, and the estimation results are as follows. During time

t \in [0, 2] s

, the position initializations of the four estimators are obtained by the trust-region method, and results as shown in Table 2. Although the optimization processing does not take clock errors and measurement errors into consideration, and hence the results are rough, it still provide a good initialization (because errors are bounded by

σ_{ρ}

) for further estimation. During time

t \in [2, 6000]

, the relative-position estimation errors and relative-velocity estimation errors are shown in Figs. 10 and 11 respectively, from which the following conclusion can be drawn. The errors of relative-position estimation are no more than 0.2 m throughout all the time after convergence and are much smaller when relative motion becomes slight, and the errors of relative-velocity estimation are no more than

0.005 m / s

. The results of the estimation for clock errors are similar, and they are omitted to save space. The estimation errors are very small, showing that the estimation is effective and helpful to controlling. Meanwhile during time

t \in [0, 6000]

, there are also four observers trying to observe

ρ_{i}

,

i = 1, \dots, 4

, and the errors of the observing results are shown in Fig. 12, from which it can be seen that the errors have been decreased to no more than 10 m before

t = 60 s

. In the control law (40)

{\hat{r}}_{i}

will be multiplied by matrix

M_{1}

and meanwhile the elements of

M_{1}

are small, then it is clear

{\hat{r}}_{i}

and

ρ_{i}

make little difference to the control law.

Fig. 8. 2D projections of motion trajectories of the four maneuvering spacecraft onto the $x - y$ plane

Fig. 9. 3D positions of the four spacecraft during maneuvering

Table 2. Optimization for Position Initialization

Relative position	Seeking origin	Optimization	True value
$ρ_{10}$	${[0, 0, 0]}^{T} m$	${[- 6.43, 734.68, - 33.28]}^{T} m$	${[- 0.04, 735.46, 0.45]}^{T} m$
$ρ_{21}$	${[0, 0, 0]}^{T} m$	${[- 11.36,735.13, - 18.65]}^{T} m$	${[- 0.11,735.46,0.45]}^{T} m$
$ρ_{32}$	${[0, 0, 0]}^{T} m$	${[26.69, 734.85, - 14.70]}^{T} m$	${[- 0.18,735.46,0.45]}^{T} m$
$ρ_{43}$	${[0, 0, 0]}^{T} m$	${[25.69,734.99,3.25]}^{T} m$	${[- 0.26,735.46,0.45]}^{T} m$

Fig. 10. Relative-position estimation errors

Fig. 11. Relative-velocity estimation errors

Conclusion

In this paper, coordinated motion problem of distributed spacecraft systems that consist of one free-flying spacecraft and several maneuvering spacecraft is considered, and a relatively complete processing framework is proposed for the coordinated motion. In the upper level, a distributed algorithm is proposed to allocate the estimators, helping offloading computations to every spacecraft. In the lower level, estimators, observers, and controllers are all designed in a distributed way. Among them, estimators are designed based on Kalman fiters and solutions of initialization problems are also presented, observers work for rough positions of maneuvering spacecraft in order to provide good approximations for controlling, and controllers work based on consensus and only use the results of estimators and observers to achieve coordinative motion. Finally, simulations are carried out where a low-Earth-orbit satellite flying scenario is taken as an instantiation. Simulation results show that estimators, observers, and controllers work well with fixed, connected, and undirected communication topologies and configuration can be formed with centimeter-level accuracy, having demonstrated the effectiveness of proposed approach.

Acknowledgments

The authors would like to thank the support by the National High Technology Research and Development Program of China (Grant No. 2013AA122601).

References

Cao, Y., and Ren, W. (2012). “Distributed coordinated tracking with reduced interaction via a variable structure approach.” IEEE Trans. Automat. Contr., 57(1), 33–48.

Crossref

Google Scholar

Cao, Y., Ren, W., and Meng, Z. (2010). “Decentralized finite-time sliding mode estimators and their applications in decentralized finite-time formation tracking.” Syst. Control Lett., 59(9), 522–529.

Crossref

Google Scholar

Chung, S. J., Ahsun, U., and Slotine, J. J. E. (2009). “Application of synchronization to formation flying spacecraft: Lagrangian approach.” J. Guid. Control Dyn., 32(2), 512–526.

Crossref

Google Scholar

Clohessy, W., and Wiltshire, R. (1960). “Terminal guidance system for satellite rendezvous.” J. Aerosp. Sci., 27(9), 653–658.

Crossref

Google Scholar

Coleman, T. F., and Li, Y. (1996). “An interior trust region approach for nonlinear minimization subject to bounds.” SIAM J. Optim., 6(2), 418–445.

Crossref

Google Scholar

Conn, A. R., Gould, N. I., and Toint, P. L. (2000). Trust region methods, Vol. 1, SIAM, Philadelphia.

Crossref

Google Scholar

Gunnam, K. K., Hughes, D. C., Junkins, J. L., and Kehtarnavaz, N. (2002). “A vision-based DSP embedded navigation sensor.” Sens. J. IEEE, 2(5), 428–442.

Crossref

Google Scholar

Horn, R. A., and Johnson, C. R. (2012). Matrix analysis, Cambridge University Press, New York.

Crossref

Google Scholar

Jadbabaie, A., Lin, J., and Morse, A. S. (2003). “Coordination of groups of mobile autonomous agents using nearest neighbor rules.” IEEE Trans. Automat. Contr., 48(6), 988–1001.

Crossref

Google Scholar

Julier, S. J., and Uhlmann, J. K. (2002). “Reduced sigma point filters for the propagation of means and covariances through nonlinear transformations.” Proc., American Control Conf., IEEE, 887–892.

Google Scholar

Lau, K., Lichten, S., Young, L., and Haines, B. (1996). “An innovative deep space application of gps technology for formation flying spacecraft.” Proc., AIAA GN&C Conf., AIAA, Reston, VA, 1–9.

Google Scholar

Lee, D., Bang, H., Butcher, E. A., and Sanyal, A. K. (2014). “Kinematically coupled relative spacecraft motion control using the state-dependent riccati equation method.” J. Aerosp. Eng., 04014099.

Crossref

Google Scholar

Li, Z., Duan, Z., and Chen, G. (2010). “Consensus of multiagent systems and synchronization of complex networks: A unified viewpoint.” IEEE Trans. Circuits Syst., 57(1), 213–224.

Crossref

Google Scholar

Meng, Z., Ren, W., and You, Z. (2010). “Distributed finite-time attitude containment control for multiple rigid bodies.” Automatica, 46(12), 2092–2099.

Crossref

Google Scholar

Meng, Z., Zhao, Z., and Lin, Z. (2013). “On global leader-following consensus of identical linear dynamic systems subject to actuator saturation.” Syst. Control Lett., 62(2), 132–142.

Crossref

Google Scholar

Mitchell, M. L. (2004). “CDGPS-based relative navigation for multiple spacecraft.” M.S. thesis, Univ. of Texas at Austin, Austin, TX.

Google Scholar

Ren, W. (2007). “Formation keeping and attitude alignment for multiple spacecraft through local interactions.” J. Guid. Control Dyn., 30(2), 633–638.

Crossref

Google Scholar

Ren, W. (2008). “On consensus algorithms for double-integrator dynamics.” IEEE Trans. Automat. Contr., 53(6), 1503–1509.

Crossref

Google Scholar

Ren, W., Beard, R. W., and Atkins, E. M. (2007). “Information consensus in multivehicle cooperative control: Collective group behavior through local interaction.” IEEE Contr. Syst. Mag., 27(2), 71–82.

Crossref

Google Scholar

Yang, T., Meng, Z., Dimarogonas, D. V., and Johansson, K. H. (2014). “Global consensus for discrete-time multi-agent systems with input saturation constraints.” Automatica, 50(2), 499–506.

Crossref

Google Scholar

Information & Authors

Information

Published In

Journal of Aerospace Engineering

Volume 29 • Issue 3 • May 2016

Copyright

This work is made available under the terms of the Creative Commons Attribution 4.0 International license, http://creativecommons.org/licenses/by/4.0/.

History

Received: Sep 26, 2013

Accepted: Jul 29, 2015

Published online: Oct 19, 2015

Discussion open until: Mar 19, 2016

Published in print: May 1, 2016

Authors

Affiliations

Xiaochu Wang, Ph.D. [email protected]

Assistant Research Scientist, Qian Xuesen Laboratory of Space Technology, China Academy of Space Technology, No. 104 Youyi Rd., Beijing 100094, China; formerly, State Key Laboratory of Precision Measurement Technology and Instruments, Dept. of Precision Instrument, Tsinghua Univ., Beijing 100084, China (corresponding author). E-mail: [email protected]; [email protected]

View all articles by this author

Kaichun Zhao, Ph.D. [email protected]

Assistant Professor, State Key Laboratory of Precision Measurement Technology and Instruments, Dept. of Precision Instrument, Tsinghua Univ., Beijing 100084, China. E-mail: [email protected]

View all articles by this author

Zheng You, Ph.D. [email protected]

Professor, State Key Laboratory of Precision Measurement Technology and Instruments, Dept. of Precision Instrument, Tsinghua Univ., Beijing 100084, China. E-mail: [email protected]

View all articles by this author

Metrics & Citations

Metrics

Citations

Download citation

If you have the appropriate software installed, you can download article citation data to the citation manager of your choice. Simply select your manager software from the list below and click Download.

Coordinated Motion Control of Distributed Spacecraft with Relative State Estimation

Abstract

Introduction

Backgrounds

Algebraic Graph Theory

Coordinates

AFF Sensor

Relative Estimation

Measurement Equations

Relative Motion Equations

Estimation

Distributed Control Based on Estimations

Analyses on Design of Control Law

Observing $ρ_{i}$ Based on $ρ_{i j}$

Design of Control Law

Simulation

Conclusion

Acknowledgments

References

Information & Authors

Information

Published In

Copyright

History

Authors

Affiliations

Metrics & Citations

Metrics

Citations

Download citation

Cited by

View Options

Media

Figures

Other

Tables

PREVIOUS ARTICLE

NEXT ARTICLE

Verify Phone

Congrats!

Abstract

Introduction

Backgrounds

Algebraic Graph Theory

Coordinates

AFF Sensor

Relative Estimation

Measurement Equations

Relative Motion Equations

Estimation

Distributed Control Based on Estimations

Analyses on Design of Control Law

Observing ρi Based on ρij

Design of Control Law

Simulation

Conclusion

Acknowledgments

References

Information

Published In

Copyright

History

Authors

Affiliations

Metrics

Citations

Download citation

Cited by

Figures

Other

Share

Copy the content Link

Share with email

Share

Request Username

Create a new account

Change Password

Password Changed Successfully

Verify Phone

Congrats!

Observing $ρ_{i}$ Based on $ρ_{i j}$