インパルス制御による円軌道上のフォーメーション
2010SE261 若園直志 2010SE262 鷲野友紀
指導教員：市川朗

1x = y = z = 0 で線形化すると,はじめにẍ − 2nẏ − 3n2 x = uxÿ + 2nẋ = uy

複数の人工衛星によるフォーメーションフライトにより, 高精度な観測や, 種々のミッションが可能となる. 飛行中の人工衛星には, 燃料が補給できないため, 低燃費の運z̈ + n z = uz用が必要となる. 2007 年に日本で打ち上げられたかぐやが得られる. (2) は HCW 方程式と呼ばれる [1],[2].

(SELENE) がフォーメーションフライトの例として挙げられる. かぐやは主衛星としておきな (OKINA)・おうな(OUNA) と呼ばれる従衛星と共にフォーメーションフライトを行い, 月の地表や, 元素・鉱物組成の調査を行った.
本研究では, フォーメーションフライトで一般的に用いられているパルス制御ではなく, インパルス制御を用いる.
パルス制御とは, 一定時間に入力を与える離散時間型の制御方法で, インパルス制御とは, 瞬間的に入力を与える離散時間型の制御方法である. 研究課題として, 燃料消費と比例関係にある総速度変化 (∆V ) を抑えた制御を行うことにより, 燃費の向上を目指す. 連続時間型のフィードバック制御を比較対象とし, インパルス制御と燃料消費の観点から比較を行う.

HCW 方程式の状態方程式は,2相対問題の方程式
半径 R0 の円軌道上の主衛星とその近傍の従衛星の相対運動を考えるため, 主衛星の重心を原点とする図 1 の回転座標系 o - {i,j,k} を考える.
従衛星
yRji

ここで

0
 0

 3n2A= 0
 0001000010002n00 −2n 000000000 −n2 000  0 10 ,B =  00 01 00とすると面内運動の解は,x(t) = 4x0 + 2y˙0 /n − (3x0 + 2y˙0 /n) cos nt+(x˙0 /n) sin nty(t) = y0 − 2x˙0 /n + (2x˙0 /n) cos nt+(6x0 + 4y˙0 /n) sin nt − (6nx0 + 3y˙0 )tẋ(t) = x˙0 cos nt + (3nx0 + 2y˙0 ) sin ntẏ(t) = (6nx0 + 4y˙0 ) cos nt−2x˙0 sin nt − (6nx0 + 3y˙0 )(4)(5)

である. (4), (5) の式をパラメータ表現すると

x(t) = 2c + a cos(nt + α)y(t) = d − 3nct − 2a sin(nt + α)z(t) = b cos(nt + β)
図 1 円軌道上の主衛星と, 運動方程式より,µ(R0 + x) + ux
R3µy + uy
R3(6)となる. ここで,このとき相対位置ベクトルを r = xi + yj + zk とするµz̈ = − 3 z + uz
R0 0 0 
0 
01

を u=[0 0 0]T とし, 面内の初期値を x0in =[x0 y0 x˙0 y˙0 ]Tz(t) = z0 cos nt + (z˙0 /n) sin ntż(t) = −nz0 sin nt + z˙0 cos nt

地球の中心ÿ = −2nẋ + n2 y −000100である. 軌道面内と軌道面外の解は独立している. 加速度R0ẍ = 2nẏ + n2 (R0 + x) −(3)外運動の解は,主衛星

oẋ = Ax + Bu, x(0) = x0

で与えられる. 面外の初期値を x0out =[z0 z˙0 ]T とする面

xr(2)

2(1)

11

a = [(3x0 + 2y˙0 /n)2 + (x˙0 /n)2 ] 2 , b = [z02 + (z˙0 /n)2 ] 2c = 2x0 + y˙0 /n, d = y0 − 2x˙0 /nsin α = −x˙0 /na,sin β = −z˙0 /nb
cos α = −(3x0 + 2y˙0 /n)/a, cos β = z0 /bである.
面内運動は c = 2x0 + y˙0 /n = 0 を満たすとき周期解となり, この条件を CW 条件と呼ぶ. フォーメーションフライが得られる. ここで, u= [ux uy uz ]1

Tは従衛星に働く加速度, R = [(R0 + x)2 + y 2 + z 2 ] 2 である. (1) を原点

トでは従衛星をこの周期解にのせる.

を用いて, uk = −Kxk となる [3]. ここで K は,インパルス制御

3

本章ではインパルスについて説明する. 状態方程式,

K = (R + B T XB)−1 B T XAdẋ = Ax + Bu, x(0) = x0

目標軌道の自由運動, 目標軌道の初期値を xf 0 とすると,を解くと,

xf k+1 = Ad xf k ,∫ tx(t) = eAt x0 +eA(t−r) Bu(r)dr(7)0となる. (7) を用いて, 時間 τ においてインパルス入力 uτは,(13)となり, (11) を用いて (13) への移行を考える. 誤差を

ek = xk − xf k とすると,
ek+1 = Ad ek − Buk ,を入れたときの応答を求めるために, パルス応答を考える.
区間 [τ, τ + h] において一定入力 uhτ を与えたときの応答

xf (0) = xf 0x(0) = x0(14)となる. ここで uk = −Kek を用いることにより, 誤差 ekは 0 に収束し, 目標の軌跡に追従する.

t < τ のとき,

本 研 究 で は (1), (2) を 軌 道 半 径 R0 お よ び 円√軌 道 のt > τ + h のとき,x(t) = eAt x0 + eAt

1h∫ τ +h

シミュレーションの準備

5x(t) = eAt x0µ/R0 3である. また時間 t を σ = t/(1/n) とおき, (x, y, z) =(x/R0 , y/R0 , z/R0 ), (ux , uy , uz ) = (1/R0 n2 )(ux , uy , uz )とすると, 状態方程式は角速度 n を用いて無次元化する. ここで n =e−Ar Buτ dr(8)τとなる. 区間 [τ, τ + h] での入力の積分は一定値 uτ とな

x′ = Ax + Bu,る. ここで h を 0 に近づけたときの極限を, 時間 τ におけx(0) = x0

るインパルス uτ の応答と呼ぶ. (8) より, インパルス入力

 0 0
 0 0
 3 0
A=
 0 0

0 00 0を加えたときの応答は,x(t) = eAt x0 + eAt e−Aτ Buτ
= eAt x0 + eA(t−τ ) Buτ(9)である. 次に一定時間ごとにインパルスを入力する. 時間kτ にインパルス入力 uk−1 を加えるとし, xk = x(kτ +)とおく. (9) を用いて,00000−10 0 0 0 
102πxk+1 = Ad xk + Buk

xk+1 = eAτ xk + Buk(10)となる. 本研究では, 1 周期に λ 回のインパルス入力を考Tえる. 周期を T とおくと τ = Tλ となる. Ad = eA λ とおくと (10) は,(11)であり, 入力を加えた直後の状態を表す漸化式が得られる.

4012000と表すことができ, (11) より, Ad = eA λ であり,

x1 = x(τ +) = eAτ x0 + Bu0xk+1 = Ad xk + Buk

100−200最適レギュレータとなる.

6シミュレーション

6.1

総速度変化

本研究では, 従衛星の初期値を xs0 =[0.01 0 0 -0.02 0.01
T0] , 目標軌道の初期値を xf 0 =[0.005 0 0 -0.01 0.005 0]Tと設定した. 本節では ∆V の変化を見るため, インパルス
入力回数を 1 周期に 4 回とする.
最適レギュレータによって設計された制御系を評価する方最適レギュレータでは, 以下のシステムと評価関数を考える.

法として, 燃料消費と比例関係にある ∆V と収束するまでの時間, 整定時間 ST を用いる. 本研究では, 整定時間は,

xk+1 = Ad xk + Buk ,J(u; x0 ) =∞∑x(0) = x0

制御軌道と目標軌道の誤差が 10−4 以下になる時間とする.
総速度変化 ∆V を面内運動 ∆Vin , 面外運動 ∆Vout に分けて表すと,TT[x (k)Qx(k) + u (k)Ru(k)] (12)k=0∆Vin =uink

k=1(12) を最小にする u(k) は, 離散時間リッカチ代数方程式,
X = Q + ATd XAd − ATd XB(R + B T XB)−1 B T XAdST∑∆Vout =ST∑
k=1|uzk |uin =√u2x + u2y表 1 入力回数の変化による ∆Vとなる. ∆V の変化をみるため, 今回最適レギュレータで用いた重み行列 Q を Q = 10I6 で固定し, 重み行列 R をR = 10r I3 の r を変化させる. このときの ∆V のグラフを
図 2 に表す.

回数

40.019∆V

50.02160.02770.03580.04290.052-3x 107.5∆Vin∆Vout

1,2 回の入力は可制御ではなくなる. 3 回の場合は収束す7るが, 目標軌道の内側を通過し, ∆V の値が大きくなるた

6.5め除外した. 5 回以上の入力を行っていくと徐々に ∆V が∆V6上がっていってしまうため, 数値的にも 4 回が最も優れて5.5いるということが確認できた.

54.56.3

重み行列43.5-1

前節で求めた入力回数を用いて, より燃料消費を抑える-0.500.5

11.52

ため, ∆V を抑える重み行列を探す. 入力に対する重み行列r

図 2 r の変化に伴うインパルス制御の ∆V の推移
図 2 より, r の増加に伴って操作量が制限され, ∆Vin

R を R = 10r I3 とし, r を −1 から 2 まで変化させつつ,状態に対する重み行列 Q を
 1を抑えて, 燃料の消費を抑えていることが分かる. 一方,0 0 0 0 0 0 0 0 p 0 0 0 
0 p 0 0 
0 0 1 00 0 0 p0
 0 1
 0 0
Q = 10 × 
 0 0

0 00 0∆Vout の変化は小さい. これは面外運動が r に依存せず,早く収束してしまうからだと考えられる. 次に, r の変化に伴う整定時間の推移を図 3 に表す.
30とし, 速度成分の p を変化させていく.

Q を変化させたときに, 各値が R の範囲で最も ∆V を抑えた値と整定時間のグラフを図 5 と図 6 に表す.

25ST20-3x 10
3.78153.77103.76

3.75

-0.5

00.5

11.52r

図 3 r の変化に伴う整定時間の推移

最小値 ∆V5-1

3.74

3.733.72

図 3 より, r の増加, つまり燃料の消費を抑えることによ

3.71

り, 収束する時間は長くなることが分かる.

6.2

3.7

入力回数

3.69
00.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.91p

図 5 Q の変化に伴う ∆V の推移

本章では, 最も ∆V を抑えられる入力回数を探す. 入力回数を 1 から 10 回まで変化させたグラフと値を図 4 と表1 に表す.
384.53643.534

整定時間

3∆V2.5322301.51280.50260246810入力回数

図 4 入力回数の変化による ∆V の推移

00.1

0.2

0.3

0.4

0.50.60.70.80.9p

図 6 Q の変化に伴う整定時間の推移

1100.064

図 5 と図 6 の結果より, ∆V を抑え, なおかつ整定時間0.024も短くする値は, p = 0.07 であり, r は r = 2 である. ま0.022

L1ΔV0.02あった. この重み行列を用いて, 連続時間型フィードバック0.018

制御との比較を行う.

0.016
L1

た, このときの ∆V は 3.70 × 10−3 で, 整定時間は 33.76 で0.014

7連続時間型フィードバック制御との比較

0.012

本節ではインパルス制御と, 連続時間型フィードバック0.01制御を比較する. 対等な条件下で比較するため, インパル

0.008

ス制御の時と同様, 連続時間型フィードバック制御に対し0.006
020406080100120STても, ノルムを抑えつつ, 整定時間も考慮した重み行列 Q

図 9 L1 ノルムと ∆V の比較を探す.

Q を変化させたときに, 各値が R の範囲で最も L1 ノルムを抑えた値と整定時間をとったグラフを図 7 と図 8 に表す.

表 2 インパルス制御と連続型制御の比較

-3

名称

面内面外

整定時間

インパルス

3.70 × 10−3
4.80 × 10−3

4.88 × 10−3
6.17 × 10−3

30.63
18.06

x 104.85

連続型制御

4.84

L1ノルム

4.834.82

表 2 より, 整定時間 5 周期以内で最も ∆V ,L1 ノルムを4.81抑えたものを比べると, インパルス制御の方が燃料消費を4.8抑えられることが表 2 からも分かる. 仮に整定時間を短く4.79することを目的とするならば, 連続型制御の方が優れてい4.78ると考えられる.

4.77

4.76
00.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

8

おわりにp

図 7 Q を変化させた時の L1 ノルムのグラフインパルス制御によるフォーメーションについて研究を行った. フォーメーションフライトは現在, 高精度な観測を行えるという利点から, 徐々に宇宙活動の際に用いられ始めている手法の一つである. そのため研究題材として40扱った. 連続時間型フィードバック制御, インパルス制御それぞれに燃料消費を抑えるために適した重み行列 Q を35設定した. それにより, 両制御が燃料消費を限りなく抑え整定時間30られるような条件でシミュレーションを行った結果, 燃料25消費を抑えられる制御は, インパルス制御であるということが分かった.

20参考文献151000.2

0.40.60.81p

図 8 Q を変化させた時の整定時間のグラフ
図 7, 図 8 より, L1 ノルムを抑えつつ, 整定時間を短くする p は p = 0.52 で r は r = 2 である. また, このときのL1 ノルムは 4.80 × 10−3 で, 整定時間は 18.06 であった.
以上の値を元に, インパルス制御と連続時間型フィードバック制御の性能を比較した結果, インパルス制御の方が燃料消費を抑えられることが分かった. また最後に収束周期を 5 周期以内と設定し, ∆V と L1 ノルムの比較を行った. 横軸を整定時間, 縦軸を ∆V と L1 ノルムとしたグラ
フと収束周期を 5 周期以内としたときの値を図 9 と表 2 に示す.

[1] M.Shibata and A.Ichikawa:Orbital Rendezvous andFlyaround Based on Null controllability with Vanishing Energy, Journal of Guidance, Control, andDynamics, Vol.30,No.4,pp.934-945,2007.
[2] B.Wie:Space Vehicle Dynamics and Control,AIAA
Education Series,Reston,Virginia,1998.
[3] Arnold, W.F.,III and A.J. Laub:Generalized Eigenproblem Algorithms and Software for Algebrai Riccati Equations, Proc. IEEE, 72(1984), pp.17461754.