解き方の全体像

最適制御問題の数値解法は、間接法と直接法の2つに大別されます。 違いは「最適である条件を先に数式で導くか」「先に問題を有限個の数で近似してしまうか」の順番です。 OpenGoddard が使う擬スペクトル法は直接法の一種です。

	間接法	直接法
考え方	先に最適性の条件を導き、それを解く（optimize, then discretize）	先に有限個の変数で近似し、普通の最適化問題として解く（discretize, then optimize）
解くもの	2点境界値問題（微分方程式）	非線形計画問題（NLP）
代表的な手法	ポントリャーギンの最小原理＋シューティング法	直接シューティング法、コロケーション法（擬スペクトル法を含む）

間接法：最適性の条件を先に導く

普通の最適化との対応

関数 $f(x)$ の最小点では $f'(x) = 0$ が成り立ちます。そこで、まず $f'(x) = 0$ という方程式を作り、それを解く、 というのが間接法の考え方です。制約付きの問題ではラグランジュ乗数を導入して同じことをします。

最適制御問題でこれに当たるのがポントリャーギンの最小原理です。運動方程式という制約の各時刻にラグランジュ乗数に当たる量 $\boldsymbol{\lambda}(t)$ を用意し、ハミルトニアン

$$H(\boldsymbol{x}, \boldsymbol{u}, \boldsymbol{\lambda}, t) = L(\boldsymbol{x}, \boldsymbol{u}, t) + \boldsymbol{\lambda}^\top \boldsymbol{f}(\boldsymbol{x}, \boldsymbol{u}, t)$$

を作ります。$\boldsymbol{\lambda}(t)$ を共状態（costate、随伴変数）といいます。最適解では次が成り立ちます。

$$\begin{aligned} \dot{\boldsymbol{x}} &= \frac{\partial H}{\partial \boldsymbol{\lambda}} = \boldsymbol{f} && \text{（運動方程式そのもの）}\\[2pt] \dot{\boldsymbol{\lambda}} &= -\frac{\partial H}{\partial \boldsymbol{x}} && \text{（共状態の方程式）}\\[2pt] \boldsymbol{u}^*(t) &= \operatorname*{arg\,min}_{\boldsymbol{u}} H && \text{（各時刻で } H \text{ を最小にする制御を選ぶ）} \end{aligned}$$

これに、境界条件と終端での $\boldsymbol{\lambda}$ の条件（横断条件）が加わります。

例：止まるまでの最短時間

最適制御問題とは の例（$\dot{x} = v,\ \dot{v} = a,\ \lvert a \rvert \le 2$、所要時間最小）に当てはめます。 所要時間は $\int_0^{t_f} 1\,dt$ なので $L = 1$ とすると、

$$H = 1 + \lambda_x v + \lambda_v a$$

共状態の方程式は $\dot{\lambda}_x = -\partial H/\partial x = 0$、$\dot{\lambda}_v = -\partial H/\partial v = -\lambda_x$ です。 よって $\lambda_x$ は定数、$\lambda_v$ は時間の1次関数です。$H$ を最小にする $a$ は

$$a^*(t) = \begin{cases} +2 & (\lambda_v(t) < 0) \\ -2 & (\lambda_v(t) > 0) \end{cases}$$

で、1次関数 $\lambda_v$ は符号が高々1回しか変わりません。したがって「最大加速から最大減速へ1回だけ切り替える」が最適だと、 計算する前にわかります。このように解の構造がわかるのが間接法の大きな利点です。

2点境界値問題とシューティング法

数値的に解く段階では、困難が現れます。状態 $\boldsymbol{x}$ は初期値が与えられていますが、共状態 $\boldsymbol{\lambda}$ の初期値はわかりません。 代わりに終端で条件が与えられます。このように始めと終わりに条件が分かれた微分方程式を2点境界値問題といいます。

典型的な解き方はシューティング法です。$\boldsymbol{\lambda}(t_0)$ を仮に決めて微分方程式を最後まで積分し、 終端条件からのずれを見て $\boldsymbol{\lambda}(t_0)$ を修正する、を繰り返します。大砲の角度を調整して的に当てるのに似ています。

左：シューティング法。初期値を決めて最後まで積分し、外れたら初期値を直します。右：コロケーション法。途中の点（ノード）の値をすべて未知数とし、全体として運動方程式と条件が成り立つように同時に動かします。

間接法の長所と短所

• 長所：解が高精度。最適解の構造（切り替えの回数など）が解析的にわかる。
• 短所：問題ごとにハミルトニアンと共状態の方程式を手で導く必要がある。 不等式制約（推力上限、高度下限など）が入ると、制約が効いている区間と効いていない区間の並び方を事前に仮定する必要がある。
• 短所：共状態の初期値には物理的な意味がないので、良い推定値を与えにくい。 しかも積分の途中で誤差が指数的に拡大しやすく、少しずれた推定値から収束しないことが多い。

直接法：先に有限個の数で近似する

直接法では、最適性の条件は導きません。代わりに、未知の関数 $\boldsymbol{x}(t), \boldsymbol{u}(t)$ を有限個の数で表し、 運動方程式・条件・目的関数をそれらの数の式に書き換えます。すると問題は

$$\min_{\boldsymbol{z}} F(\boldsymbol{z}) \quad \text{subject to}\quad \boldsymbol{c}(\boldsymbol{z}) = \boldsymbol{0},\quad \boldsymbol{d}(\boldsymbol{z}) \ge \boldsymbol{0}$$

という、有限個の変数 $\boldsymbol{z}$ についての制約付き最適化問題になります。これを非線形計画問題（NLP）といいます。 NLP は汎用のソルバーで解けます。OpenGoddard は SciPy の SLSQP（逐次二次計画法）を使います。 逐次二次計画法は、目的関数を2次式、制約を1次式で近似した簡単な問題を繰り返し解いて、解に近づく方法です。

直接法には、関数の表し方によっていくつかの種類があります。

直接シューティング法

制御 $\boldsymbol{u}(t)$ だけを有限個の値（たとえば区分的に一定）で表し、状態は数値積分で計算します。 変数は少なくて済みますが、間接法のシューティングと同じく、積分の途中で誤差が拡大しやすいという弱点があります。

コロケーション法

制御だけでなく状態も、時間方向の点（ノード）での値を未知数にします。運動方程式は数値積分するのではなく、 「各ノードで、状態の時間微分が運動方程式の右辺に一致する」という等式制約として課します（上の図の右側）。

最初は運動方程式を満たしていない点の並びから始め、最適化の過程で目的関数の改善と運動方程式の充足を同時に進めます。 長い時間の積分を何度も繰り返さないので、初期推定が悪くても収束しやすくなります。

コロケーション法は、ノード間を何で結ぶかでさらに分かれます。

• 局所コロケーション：短い区間ごとに低次の多項式（1次や3次）でつなぐ。台形則や Hermite–Simpson 法など。 素直で頑丈だが、精度を上げるには多くの点が必要で、点の数に対して誤差はゆっくり（べき乗で）しか減らない。
• 大域コロケーション＝擬スペクトル法：区間全体を1本の高次多項式で表し、ノードを特別な位置に置く。 解が滑らかなら、点の数を増やすと誤差が指数的に減る。

直接法の長所と短所

• 長所：ハミルトニアンや共状態を導く必要がなく、運動方程式・条件・目的関数を書けば解ける。
• 長所：不等式制約をそのまま書ける。制約が効く区間を事前に仮定しなくてよい。
• 長所：初期推定は状態と制御という物理的な量で与えればよく、見当をつけやすい。
• 短所：変数の数が多い（ノード数 × 状態・制御の数）。
• 短所：得られるのは離散化した問題の解で、精度は離散化の細かさに依存する。また局所的な最適解である。

擬スペクトル法の立ち位置

擬スペクトル法は「直接法」の中の「大域コロケーション法」です。主な利点は次の3つです。

1. 少ない点で高い精度。 滑らかな解なら、数十点程度で非常に高い精度が得られます。 例題03（最速降下）では、ノード数を 6、8、10 と増やすと所要時間の誤差が およそ $10^{-4}$、$10^{-7}$、$10^{-11}$ と急速に小さくなり、12 ノードで計算機の丸め誤差の程度に達します（擬スペクトル法）。
2. 問題が小さく、速く解ける。 点が少ないので NLP の変数も少なく、小型の汎用ソルバー（SLSQP）でも扱えます。
3. 間接法の情報も得られる。 NLP のラグランジュ乗数から、間接法の共状態 $\boldsymbol{\lambda}(t)$ を推定できることが知られています （Fahroo & Ross, 2001）。直接法で解いた結果を、最小原理で検算できるということです。

一方、不連続を含む解は苦手です。1本の多項式で、急に切り替わる推力や段分離での質量の跳びを表そうとすると、 切り替えの前後で振動が生じます。対策は、時間を複数の区間に分けて、区間ごとに別の多項式を使うことです （例題07 と 例題08 の比較を参照）。

擬スペクトル法にはノードの選び方によって LG 法（Gauss）、LGR 法（Radau）、LGL 法（Lobatto）などの種類があります。 OpenGoddard は区間の両端をノードに含む LGL 法を使います。両端がノードなので、境界条件や区間の接続を素直に書けます。 次のページで、その仕組みを順を追って説明します。