JP6991326B2 - 線形計画問題求解システム、解候補算出装置、最適解算出装置、宇宙機のスラスタ制御装置および飛翔体制御装置並びに線形計画問題の求解方法 - Google Patents

線形計画問題求解システム、解候補算出装置、最適解算出装置、宇宙機のスラスタ制御装置および飛翔体制御装置並びに線形計画問題の求解方法 Download PDF

Info

Publication number
JP6991326B2
JP6991326B2 JP2020520941A JP2020520941A JP6991326B2 JP 6991326 B2 JP6991326 B2 JP 6991326B2 JP 2020520941 A JP2020520941 A JP 2020520941A JP 2020520941 A JP2020520941 A JP 2020520941A JP 6991326 B2 JP6991326 B2 JP 6991326B2
Authority
JP
Japan
Prior art keywords
constraint
vector
solution
dual
optimum
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2020520941A
Other languages
English (en)
Other versions
JPWO2019224954A1 (ja
Inventor
雄太 浅野
岳宏 西山
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Mitsubishi Electric Corp
Original Assignee
Mitsubishi Electric Corp
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Mitsubishi Electric Corp filed Critical Mitsubishi Electric Corp
Publication of JPWO2019224954A1 publication Critical patent/JPWO2019224954A1/ja
Application granted granted Critical
Publication of JP6991326B2 publication Critical patent/JP6991326B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B64—AIRCRAFT; AVIATION; COSMONAUTICS
    • B64G—COSMONAUTICS; VEHICLES OR EQUIPMENT THEREFOR
    • B64G1/00—Cosmonautic vehicles
    • B64G1/22—Parts of, or equipment specially adapted for fitting in or to, cosmonautic vehicles
    • B64G1/24—Guiding or controlling apparatus, e.g. for attitude control
    • B64G1/244—Spacecraft control systems
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B64—AIRCRAFT; AVIATION; COSMONAUTICS
    • B64G—COSMONAUTICS; VEHICLES OR EQUIPMENT THEREFOR
    • B64G1/00—Cosmonautic vehicles
    • B64G1/22—Parts of, or equipment specially adapted for fitting in or to, cosmonautic vehicles
    • B64G1/24—Guiding or controlling apparatus, e.g. for attitude control
    • B64G1/244—Spacecraft control systems
    • B64G1/245—Attitude control algorithms for spacecraft attitude control
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B64—AIRCRAFT; AVIATION; COSMONAUTICS
    • B64G—COSMONAUTICS; VEHICLES OR EQUIPMENT THEREFOR
    • B64G1/00—Cosmonautic vehicles
    • B64G1/22—Parts of, or equipment specially adapted for fitting in or to, cosmonautic vehicles
    • B64G1/24—Guiding or controlling apparatus, e.g. for attitude control
    • B64G1/26—Guiding or controlling apparatus, e.g. for attitude control using jets
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B64—AIRCRAFT; AVIATION; COSMONAUTICS
    • B64U—UNMANNED AERIAL VEHICLES [UAV]; EQUIPMENT THEREFOR
    • B64U50/00—Propulsion; Power supply
    • B64U50/10—Propulsion
    • B64U50/18—Thrust vectoring
    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05D—SYSTEMS FOR CONTROLLING OR REGULATING NON-ELECTRIC VARIABLES
    • G05D1/00—Control of position, course, altitude or attitude of land, water, air or space vehicles, e.g. using automatic pilots
    • G05D1/08—Control of attitude, i.e. control of roll, pitch, or yaw
    • G05D1/0808—Control of attitude, i.e. control of roll, pitch, or yaw specially adapted for aircraft
    • G05D1/0858—Control of attitude, i.e. control of roll, pitch, or yaw specially adapted for aircraft specially adapted for vertical take-off of aircraft
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/12—Simultaneous equations, e.g. systems of linear equations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B64—AIRCRAFT; AVIATION; COSMONAUTICS
    • B64U—UNMANNED AERIAL VEHICLES [UAV]; EQUIPMENT THEREFOR
    • B64U2201/00—UAVs characterised by their flight controls

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Remote Sensing (AREA)
  • Chemical & Material Sciences (AREA)
  • Combustion & Propulsion (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Automation & Control Theory (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Operations Research (AREA)
  • Computing Systems (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Complex Calculations (AREA)

Description

この発明は、ある入力ベクトルに対して、線形制約を満たしつつ線形な評価関数を最小にする出力ベクトルを算出する線形計画問題求解システムおよびこれを用いた制御装置に関するものである。
線形計画問題の代表的な解法として、内点法とシンプレックス法がある。内点法は，反復計算が含まれるため、計算負荷が大きく、収束に至るまでの時間が見積もりにくい。また、シンプレックス法は、有限回の反復で計算が終了するものの、収束に至るまでの時間が見積もりにくいというデメリットが存在する。
特許文献1は、上限下限制約を持つ線形計画問題に対して、変数の定義域の上限値と変数との差に変数を置き換えることによって、上限と下限とを逆転させ、内点法の限界を回避して解の計算を続行することを開示する。
また、上記線形計画問題の解法を用いる対象として、複数のアクチュエータを制御するものが考えられる。例えば、複数のスラスタの噴射によって姿勢制御する宇宙機のスラスタ制御、ロータの推力を一方向にしか出力できないドローンなどである。宇宙機のスラスタの噴射量は、負の値を採らず、噴射燃料最小化は、線形計画問題に帰着する。ドローンのロータ制御においては、推力の大きさよりもトルクを優先して姿勢を迅速に変更する場合、ロータ回転数の上下限制約のもとで，可能な限り目標となるトルクを実現する問題が、線形計画問題となる。
特許文献２は、スラスタ駆動において、予め計算し保管したデータを基に、最小燃料の駆動ベクトルを求める線形計画問題をシンプレックス法で高速に解く方法が記載される。この方法で制御ループ内に最適解が得られなかった場合には，最適性若しくは制約条件を無視して解を出力する、またはスラスタ駆動をしないとする方法を開示する。
特開２００１－１８４３３４ 特表２００２－５０５６３２
従来の線形計画問題求解システムでは、反復を伴う内点法やシンプレックス法を用いるため、最適解が得るまでに必要な時間の見積もりが困難で、制御サイクル内で最適解が得られない可能性がある。制御サイクル内に最適解が得られない場合は、最適性が犠牲になるまたは制約条件を満たさなくなる問題があった。
本発明の線形計画問題求解システムは、与えられた制約係数行列、随時変化する制約ベクトルおよび制約ベクトルの次元数である制約次元数と異なる変数次元数の変数を持つ変数ベクトルの連立式で表現される制約条件の下で変数ベクトルおよび費用ベクトルの関数である目的関数を最適化する線形計画問題の解を求める線形計画問題求解システムにおいて、制約係数行列および費用ベクトルの情報を入力として、線形計画問題を主問題とする双対問題および双対問題の制約条件の式のうち等号が成立する組合せを表す有効セットを全て求め、有効セットごとに制約条件が有効となる実行可能な双対解候補を求めて有効セットと関連付けて双対解候補を記憶部に記憶する双対解候補探索部と、制約ベクトルを入力とし、記憶部に記憶された双対解候補と制約ベクトルとの内積から最適となる有効セットを選択し、選択した有効セットに対応する実行可能基底解を最適解として求めて出力する最適解算出装置とを備えるものである。
この発明は、最適性を犠牲にする、または制約条件を満たさないことなく、定められた時間内で確実に線形計画問題の解を求めることができる。
本発明の実施の形態１に係る線形計画問題求解システムの構成を説明する図である。 本発明の実施の形態１に係る解候補算出装置の構成を説明する図である。 本発明の実施の形態１に係る最適解算出装置の構成を説明する図である。 本発明の実施の形態２に係る解候補算出装置の構成を説明する図である。 本発明の実施の形態２に係る最適解算出装置の構成を説明する図である。 本発明の実施の形態４に係る宇宙機のスラスタ制御装置の構成を説明する図である。 本発明の実施の形態５に係る解候補算出装置の構成を説明する図である。 本発明の実施の形態５に係る最適解算出装置の構成を説明する図である。 本発明の実施の形態６に係る解候補算出装置の構成を説明する図である。 本発明の実施の形態６に係る最適解算出装置の構成を説明する図である。 本発明の実施の形態７に係る飛翔体のロータ制御装置の構成を説明する図である。 本発明の実施の形態に係る解候補算出装置のハードウェア構成を説明する図である。 本発明の実施の形態に係る最適解算出装置のハードウェア構成を説明する図である。
図１は、この発明の実施の形態１による線形計画問題求解システムの構成を示す。線形計画問題求解システム１は、入力に対して、線形制約を満たしつつ線形な評価関数を最小または最大にする出力ベクトルを算出するものである。図において、線形計画問題求解システム１は、対象とする線形計画問題を繰り返し解く際に、不変となる情報および随時変化する変化情報を入力として、入力された線形計画問題の最適解を求めて、最適解を出力するシステムである。線形計画問題求解システム１は、解候補算出装置２及び最適解算出装置３により構成される。
線形計画問題求解システム１は、与えられた制約係数行列、随時変化する制約ベクトルおよび制約ベクトルの次元数である制約次元数と異なる変数次元数の変数を持つ変数ベクトルの連立式で表現される制約条件の下で変数ベクトルおよび費用ベクトルの関数である目的関数を最適化する線形計画問題の解を求める線形計画問題求解システムである。
線形計画問題求解システム１は、制約係数行列および費用ベクトルの情報を入力として、上記線形計画問題を主問題とする双対問題およびこの双対問題の制約条件の式のうち有効な組合せを表す有効セットを全て求め、求めた有効セットごとに制約条件が有効となる実行可能な双対解候補を求めて有効セットと関連付けて双対解候補を記憶部に記憶する双対解候補探索部と、制約ベクトルを入力とし、記憶部に記憶された双対解候補と制約ベクトルとの内積から最適となる有効セットを選択し、選択した有効セットに対応する実行可能基底解を最適解として求めて出力する最適解算出装置とを備える。
さらに、双対解候補探索部は、双対問題が変数次元数個の式を持つ双対制約条件を含む問題として求め、最適解算出装置は、制約係数行列のうち最適有効セットの情報を使って求めた行列の逆行列と制約ベクトルとの積から変数ベクトルの最適解を求める最適解演算部を備える。ここで、有効セットは、双対制約条件の複数の式のうち変数次元数および制約次元数の少ない方の次元数個の有効な組合せを式の識別子で表したものであり、双対解候補は、有効セットごとに識別子が表す前記双対制約条件の式を満たす解ベクトルである。
解候補算出装置２は、対象とする線形計画問題の不変となる情報を入力として、入力された情報から対象とする線形計画問題を求め、求めた線形計画問題を主問題とする双対問題を求め、求めた双対問題の不等式制約条件の一部が有効となる実行可能な解候補と解候補に対応する有効セットをすべて求めて記憶部に記憶させる。 ここで「不等式制約条件の一部が有効」とは、主問題の制約条件の数と変数ベクトルの次元数の内の少ない方の数となる。
最適解算出装置３は、随時変化する変化情報を入力とし、解候補算出装置２で求めた解候補と入力の内積から最適となる有効セットを選択し、有効セットに対応する実行可能基底解を最適解として求めて出力する。
解候補算出装置２は、上記不変となる情報を入力する問題入力部２１、および問題入力部２１で入力された情報から線形計画問題の双対問題を求め、求めた双対問題の双対制約条件の内の有効な組合せである有効セットの内の解が一意に求められる有効セットを全て求め、求めた有効セットごとに双対制約条件を満たす解を双対解候補として全て求めて有効セットと関連付けて双対解候補を記憶部に記憶する双対解候補探索部２２によって構成される。解候補算出装置２は、有効セット双対解記憶部４に格納するための情報を算出する。
双対解候補探索部２２は、入力された情報から、線形計画問題の双対問題を求める双対問題作成部２２１、求めた双対問題の双対制約条件の内の有効な組合せである有効セットの内の解が一意に求められる有効セットを全て求める全有効セット作成部２２２、および求めた有効セットごとに双対制約条件を満たす解を双対解候補として全て求めて有効セットと関連付けて双対解候補を記憶部に記憶する双対解算出部２２３によって構成される。
最適解算出装置３は、線形計画問題の随時変化する変化情報を入力する制御ベクトル入力部３１、制御ベクトル入力３１で入力された変化情報から、解候補探索部２２によって求めた記憶部に記憶された双対解候補から双対問題が最適となる有効セットである最適有効セットを求める最適有効セット演算部３２、求めた最適有効セットと不変となる情報と変化情報とから線形計画問題の最適解を求めて出力する最適解演算部３３により構成される。
ここで、線形計画問題求解システム１が最適化する線形計画問題として、非退化な標準形の線形計画問題を例にとり説明する。ここで扱う線形計画問題は、与えられた制約係数行列A、随時変化する制約ベクトルｂ及び制約ベクトルｂの次元数である制約次元数Mより多い変数次元数Ｎの変数を持つ非負の変数ベクトルｘの連立方程式で表される制約条件の下で、変数ベクトルｘの関数である目的関数を最適化（最小に）する問題で、以下のいわゆる標準形の線形計画問題の式（１）から（３）で表現される。
Figure 0006991326000001
 ・・・・・・（１）
Figure 0006991326000002
 ・・・・・・（２）
Figure 0006991326000003
 ・・・・・・（３）
ここで、費用ベクトルｃは、変数次元数Nの次元数を持つN次元（列）ベクトル、制約係数行列Aは、制約次元数M×変数次元数NのM×N次元行列、制約ベクトルｂは、制約次元数MのM次元（列）ベクトルである。また、変数ベクトルｘは変数次元数NのN次元（列）ベクトルであり、変数次元数Nは、制約次元数Mより大きいものとする。また、Ｔは、転置行列を表す。ベクトルはすべて列ベクトルを表し、費用ベクトルｃの転置行列cTは、行ベクトルとなる。式（３）のs.t.は、such thatを表し、つづく条件を満たすことを表す。
式（２）は、制約係数行列Ａと変数ベクトルｘの積が制約ベクトルｂとなること、式（３）は、変数ベクトルｘのすべての要素が０以上であることを表す。
式（１）は、式（２），（３）を満たす変数ベクトルｘのうちで，費用ベクトルｃと変数ベクトルｘの内積を最小にする変数ベクトルｘを求めることを意味する．
一般的な線形計画問題は、上述の制約係数行列Ａが与えられ、随時変化する制約ベクトルｂの制約次元数Ｍよりも変数ベクトルｘの変数次元数Ｎが大きい問題において、制約係数行列Ａと変数ベクトルｘとの積が制約ベクトルｂとなり、変数ベクトルｘが非負である制約のもとで、与えられる任意の実数を要素とする費用ベクトルｃと変数ベクトルｘとの内積を最小にする問題を対象とする。ただし、制約係数行列Ａおよび費用ベクトルｃが不変である。
上記問題について、例えば、制御対象であるアクチュエータの制御変数を最適化する問題に適用すると、制約係数行列Ａは、アクチュエータの幾何学的配置等に依存する変換行列となる。変数ベクトルｘは、アクチュエータの制御変数である。また、費用ベクトルｃは、制御変数ｘの重みとなる。ここで、特徴は、変換行列である制御係数行列Ａと制御変数である変数ベクトルｘの重みである費用ベクトルｃとが不変な情報として与えられ、逐次変化する力やトルクなどの目標値となる制御ベクトルｂに応じて、力やトルクなどの目標値のパラメータ数よりも多い数の制御変数である変数ベクトルｘを最適化する問題を対象とすることである。
さらに、宇宙機のスラスタを制御対象とする場合には、制約係数行列Ａは、スラスタの幾何学的配置のみに依存する行列、費用ベクトルｃは、すべての要素が１のベクトル、制約ベクトルbは、目標並進力およびトルクを並べたベクトル、変数ベクトルｘは、各スラスタの噴射量となる。さらに、並進３自由度および回転３自由度を制御する場合には、制約次元数Mは６となり、変数次元数Nはスラスタの本数を表すことになる。
なお、あらゆる線形計画問題は、式（１）～（３）に示す標準形の線形計画問題に変換して帰着できるから、他の形で表現される線形計画問題を対象とする場合であっても、変換して上記標準形の線形計画問題にして最適化しても良い。また、退化した問題である場合も摂動を加えることにより非退化な問題に変形することができるため、退化している問題にも適用できる（退化、非退化の定義については双対問題の部分で後述する）。さらに、ここで扱う問題は、式（１）～（３）で示される標準形の線形計画問題の制約係数行列Aと費用ベクトルｃは、変化せず、制約ベクトルｂが時間的に変化する問題を繰り返し解く場合に、高速に求解を行う方法について具体的に説明する。
問題入力部２１は、上記制約係数行列Ａ、費用ベクトルcを入力する。具体的には、例えば、表形式のユーザインタフェースでユーザが入力するようにしても良いし、制約係数行列Ａ、費用ベクトルcの各要素の情報を持つファイルを入力するようにしても良い。
双対問題作成部２２１は、問題入力部２１で入力された制約係数行列Ａ、費用ベクトルcから線形計画問題を作成し、この線形計画問題を主問題とする双対問題を作成する。
双対問題作成部２２１が作成する双対問題は、変数次元数Ｎ及び制約次元数Ｍのうち大きい方の次元数個の不等式を持つ双対制約条件を含む情報として求められる。例えば、上記非退化の式（１）～（３）で示される標準形の線形計画問題を主問題としたときの双対問題は、下記のように表現できる．
Figure 0006991326000004
 ・・・・・・（４）
Figure 0006991326000005
 ・・・・・・（５）
ここで、不等式制約の式（５）は、制約係数行列Aの転置行列と双対問題の解ベクトルｙとの積が、費用ベクトルｃ以下であることを意味する。式（４）は、不等式制約の式（５）を満たす解ベクトルｙのうちで、制約ベクトルｂと解ベクトルｙの内積を最大にするものを求めることを意味する。
上記の双対問題において、制約係数行列Aと費用ベクトルcとが固定であるという上述の問題設定から、式（５）の不等式制約をあらかじめ解くことができる。
不等式制約式（５）を満足する解ベクトルｙの領域である実行可能領域はＭ次実数空間内の凸多面体となり，線形計画問題の一般的な性質から、不等式制約式（５）を満足する解ベクトルｙの実行可能領域の頂点の一つが必ず最適解となる。この頂点は、実行可能基底解と呼ばれている。
図１２は、本実施の形態の解候補算出装置２のハードウェア構成の例を示す図である。図において、入出力インターフェース１００５は、問題入力部２１の入力、ユーザインタフェースのための表示、並びに求めた有効セット、双対解候補および有効セットと双対解候補との関係を最適解算出装置３または外部記憶装置への出力を行い、プロセッサ１００２は、双対解候補探索部２２が行う線形計画問題の双対問題、有効セット、および双対解候補を求める演算を行い、メモリ１００３は、解候補算出装置２で実行するプログラムを保存し、双対解候補探索部２２で求めた双対解候補と有効セットとを関連付けて記憶して、有効セット双対解記憶部４として機能する。また、バス１００４は、プロセッサ１００２、メモリ１００３、入出力インターフェース１００５の間の情報、および信号を伝達する。なお、有効セット双対解記憶部４は、解候補算出装置２のメモリ１００３としても良いし、外部の記憶装置であっても良いし、最適解算出装置３のメモリであっても良い。
図１３は、本実施の形態の最適解算出装置３のハードウェア構成の例を示す図である。図において、入出力インターフェース１１０５は、線形計画問題の随時変化する変化情報を入力する制御ベクトル入力部３１の入力、双対解候補探索部２２で求めた有効セットと関連付けられた双対解候補の入力を行い、プロセッサ１１０２は、最適有効セット演算部３２が行う最適有効セットを求める演算および最適解演算部３３が行う最適解を求める演算を行い、メモリ１１０３は、最適解算出装置３で実行するプログラムを保存し、双対解候補探索部２２で求めた有効セットと関連付けられた双対解候補を記憶する。また、バス１１０４は、プロセッサ１１０２、メモリ１１０３、入出力インターフェース１１０５の間の情報、および信号を伝達する。なお、有効セット双対解記憶部４は、最適解算出装置３のメモリ１１０３としても良いし、外部の記憶装置であっても良いし、解候補算出装置２のメモリ１００３であっても良い。
また、上記メモリ１００３、１１０３は、例えば、RAM、ROM、フラッシュメモリー、EPROM、EEPROM等の、不揮発性または揮発性の半導体メモリや、磁気ディスク、フレキシブルディスク、光ディスク、コンパクトディスク、ミニディスク、ＤＶＤ等が該当する。
また、図１３は、上記で例として記載した、制御対象であるアクチュエータの制御変数を最適化した問題に適用した場合の、最適解算出装置３を含む装置全体の構成の例をも表す。図において、入出力インターフェース１１０５に接続する受信装置１１０６は、逐次変化する力やトルクなどの目標値となる制御ベクトルｂを受信し、入出力インターフェース１１０５に接続するアクチュエータ制御装置１１０７は、最適解算出装置３で求めた最適化された変数ベクトルｘ、すなわちアクチュエータの制御変数によって、アクチュエータ１１０８を制御する。
解候補算出装置２は、双対問題の実行可能基底解を最適解候補としてすべて算出する。最適解となる実行可能基底解は、下記に示す手順で求めることができる。
図２に解候補算出装置２の処理のフローチャートを示す。
解候補算出装置２のＳｔｅｐ１は、双対問題を作成する。制約係数行列Ａは、制約次元数Ｍ×変数次元数ＮのM×N次元行列であるから、双対問題の不等式制約は、変数次元数N個の不等式となる。このステップは、双対問題作成部２２１で実行される。
解候補算出装置２のＳｔｅｐ２は、変数次元数Ｎ個の不等式制約のうち、制約次元数Ｍ個の制約の等号が成立する、すなわち有効になるとして、有効となる制約の組み合わせを有効セットＩとする。有効セットＩは、変数次元数Ｎ個の不等式制約のうち、有効となる制約の組み合わせである。ここで、不等式制約が有効であるとは、等号が成立する、すなわち等式となることを意味する。例えば、変数次元数Ｎ = 4として、１番目と３番目の不等式制約が有効である場合、有効セットＩは、{１，３}となる。
解候補算出装置２のＳｔｅｐ３は、制約係数行列Ａの第ｉ列の制約次元数MのM次元ベクトルをｉ列制約係数ベクトル ａｉ として、ｉ列制約係数ベクトルａｋ （ｋ∈Ｉ）をｋの小さい順に横に並べた制約次元数MのM×M行列を基底行列 ＡＩにする。さらに、基底行列 ＡＩ のランクが、制約次元数M、すなわちフルランクであるかどうかを確認する。フルランクである場合は、Ｓｔｅｐ４へ進み、フルランクでない場合は、有効セットを記憶することなく、Ｓｔｅｐ１へ戻り次の有効セットに移行する。
解候補算出装置２のＳｔｅｐ４は、下記式（６）のｙに関する連立方程式の解を双対基底解ｙｂ Ｉ として算出する。
Figure 0006991326000006
 ・・・・・・（６）
ここで基底費用ベクトルｃＩは、費用ベクトルcのｋ番目の要素（ｋ∈Ｉ）をｋの小さい順に並べた制約次元数MのM次元ベクトルである。また、基底行列AＩの転置行列を転置基底行列AＩ Tにする。具体的には、ガウスの消去法などを用いて双対基底解ｙｂ Ｉを求める。
解候補算出装置２のＳｔｅｐ５は、Ｓｔｅｐ３で得られた双対基底解ｙｂ Ｉ が、不等式制約条件 式（５）を満たしているか確認する。不等式制約条件 式（５）を満たしている場合は、Ｓｔｅｐ６に進み，満たしていない場合は、Ｓｔｅｐ１へ戻り次の有効セットに移行する．
解候補算出装置２のＳｔｅｐ６は、有効セットＩと対応する双対基底解 ｙｂ Ｉを双対解候補 ｙＩとして有効セット双対解記憶部４に保存する。さらに具体的には、有効セットＩと、Ｓｔｅｐ４で求めた双対解候補 ｙＩとを関連付けて、有効セット双対解記憶部４に記憶する。その後、解候補算出装置２のＳｔｅｐ１へ戻り、次の有効セットＩへ移行する。
有効セットＩが異なっても、言い換えると制約次元数Ｍ個の制約の等号が成立する不等式が異なっても、双対基底解ｙｂ Ｉが一致することがある。このような基底解は退化しているといわれる。また、双対解候補ｙＩの中に、退化した基底解が一つでも存在するときには、線形計画問題は退化しているといわれる。ここで、線形計画問題が、非退化であるとする。すると、有効セットＩと双対解候補ｙＩとは、一対一に対応する。
次に、最適解算出装置３について、具体的に述べる。最適解算出装置３は、変化する制約ベクトルｂが、時間ごとに繰り返し与えられて、制約ベクトルｂに応じた最適解ｘ＊を高速に算出する装置である。図３に最適解算出装置３のフローチャートを示す．
最適解算出装置３のＳｔｅｐ１は、有効セット双対解記憶部４に記憶される、すべての有効セットＩに対して、制約ベクトルｂと双対解候補ｙＩとの内積を計算する。
最適解算出装置３のＳｔｅｐ２は、Ｓｔｅｐ１で求めた内積が、最大となる有効セットＩを最適有効セットＩ＊ と定義し、以下の式（７）によって、最適有効セットＩ＊を求める。
Figure 0006991326000007
 ・・・・・・（７）
式（７）を満たす最適有効セットＩ＊が、複数存在する場合は、どれを選んでも主問題の最適解は変化しないため、その中の一つだけを任意に選んで良い。
最適解算出装置３のＳｔｅｐ３は、標準形の線形計画問題の式（２）について、最適有効セットＩ＊を基底とした実行可能基底解を最適解ｘ＊ として求める。実行可能基底解の算出には、次式の線形連立方程式を解くことになる。
Figure 0006991326000008
 ・・・・・・（８）
式（８）は、最適有効セットＩ＊について、制約次元数Ｍ×Ｍの制約係数行列Aと、制約次元数Ｍ次の列ベクトルである最適解ｘ＊ との積が、制約ベクトルｂとなることを表す主問題の制約条件式（２）を設定した連立方程式である。この連立方程式は、例えば、ガウスの消去法を用いて式（８）の連立方程式を解いて、最適解ｘ＊ を求める。
以上によって、主問題に対する制約ベクトルｂについての最適解ｘ＊を算出することができる。
本実施の形態によれば、与えられた線形計画問題を主問題とする双対問題を求め、双対問題の制約条件が時間的に変化しないことを利用して、前記制約条件が有効となる組み合わせである有効セットＩごとに、前記双対問題の最適解候補をオフラインで算出し、前記有効セットとそれに対応する前記最適解候補を有効セット双対解記憶部に記憶したため、最適な有効セットである最適有効セットＩ＊、および対応する基底解を最適解ｘ＊として高速に算出することができる。
本実施の形態によれば、制約係数行列Ａ、費用ベクトルｃが固定であり、制約ベクトルｂが変化する非退化な標準形線形計画問題求解システムにおいて、制約ベクトルｂに依存しない双対解候補ｙＩおよび前記双対解候補に対応する有効セットＩを探索し、有効セット双対解記憶部４に保存する解候補算出装置２と、有効セット双対解記憶部４に保存された双対解候補ｙＩと制約ベクトルｂの内積を計算し、最大となる有効セットＩ＊を基底とした主問題の実行可能基底解を最適解ｘ＊として算出する最適解算出装置３とを備えた構成となっている。このように構成したので、双対解候補ｙＩを予めオフラインで計算しておくことができ、時間ごとに変化する制約ベクトルｂが与えられてから、主問題の最適解ｘ＊を算出するまでの過程に反復計算が不要となる。このことによって、計算の高速化および信頼性の向上を実現することができる。また、最適解を算出するまでの時間を見積もることができるため、制御サイクル内に最適解を確実に得られ、最適性および制約条件を犠牲にすることなく実時間で制御することができる。
実施の形態２．
実施の形態２に係る線形計画問題求解システム１は、解候補算出装置２、最適解算出装置３、および有効セット双対解記憶部４から構成される点で、実施の形態１と同様であるが、解候補算出装置において有効セット双対解記憶部４に格納する情報を追加することで、さらに求解を高速化することができる。具体的には、実施の形態２では、図４に示すように、図２に示した解候補算出装置２のＳｔｅｐ５、および最適解算出装置３のＳｔｅｐ３を下記のように修正する。
図４は、解候補算出装置２の処理の流れを示す図である。図において、解候補算出装置２のＳｔｅｐ５’は、有効セットＩに対応する双対基底解 ｙｂ Ｉを双対解候補ｙＩとして、有効セット双対解記憶部４に保存する（以上は、実施の形態１のＳｔｅｐ５に同じ）。さらに、有効セットＩに対する制約次元数Ｍ×Ｍの制約係数行列Ａの逆行列である分配行列ＡＩ －１を算出し、有効セットＩごとに関連付けて、有効セット双対解記憶部４に分配行列ＡＩ －１を保存する。有効セットＩに対する制約次元数Ｍ×Ｍの制約係数行列Ａの逆行列である分配行列ＡＩ －１を求めるには、例えば、制約係数行列ＡをＬＵ分解すれば良い。その後、解候補算出装置２のＳｔｅｐ１へ戻り、次の有効セットへ処理を移行する。
実施の形態１では双対解候補ｙＩを探索して終了したが、各最適解（有効セット）に対応した主問題の分配行列ＡＩ －１もあらかじめ計算しておく。また、実施の形態１の最適解算出装置３において、図５に示すように最適解算出装置３のＳｔｅｐ３を下記のように修正する。
図５は、最適解算出装置３の処理の流れを示す図である。図において、最適解算出装置３のＳｔｅｐ３’は、有効セット双対解記憶部４に記憶されている最適有効セットＩ＊に関連付けられた、分配行列ＡＩ －１を読出し、最適有効セットの分配行列ＡＩ* －１とし、制約次元数Ｍ×Ｍの最適有効セットの分配行列と制約ベクトルｂの積ＡＩ －１ｂを求め、最適解ｘI* を計算する。Ｓｔｅｐ３’では、既に求まっている分配行列ＡＩ* －１を用いるため、逆行列の計算をする必要がなく、最適解算出装置３の演算量をさらに少なくできる。
なお、解候補算出装置２は、実施の形態１に記載の図１２と同様に、プロセッサ１００２、メモリ１００３、入出力インターフェース１００５のハードウェア構成で実現できる。また、プロセッサ１００２は、制約係数行列Ａの逆行列である分配行列ＡＩ －１を算出し、メモリ１００３は、有効セットＩごとに関連付けて、有効セット双対解記憶部４に分配行列ＡＩ －１を記憶する。
また、最適解算出装置３は、実施の形態１に記載の図１３と同様に、プロセッサ１１０２、メモリ１１０３、入出力インターフェース１１０５のハードウェア構成で実現できる。また、メモリ１１０３は、有効セットＩごとに関連付けて、有効セット双対解記憶部４に分配行列ＡＩ －１を記憶し、プロセッサ１１０２は、最適有効セットＩ＊に関連付けられた、分配行列ＡＩ －１を読出し、分配行列と制約ベクトルｂの積ＡＩ －１ｂを求め、最適解ｘI* を計算する。
本実施の形態の線形計画問題求解システム１によれば、解候補算出装置２で有効セットＩに対する制約次元数Ｍ×Ｍの制約係数行列Ａの逆行列である分配行列ＡＩ －１を算出し、有効セットＩごとに関連付けて有効セット双対解記憶部４に記憶し、最適解算出装置３が、有効セット双対解記憶部４に記憶されている最適有効セットＩ＊に関連付けられた、分配行列ＡＩ －１を読出して最適解ｘI* を計算するため、最適解算出装置３の演算量をさらに少なくできる。これによって、制御サイクルをさらに短くして、応答が早いまたは滑らかな制御を実現することができる。
実施の形態３．
上記実施の形態では、制約係数行列Ａおよび費用ベクトルｃが変化せず、固定として扱ったが、制約係数行列Ａおよび費用ベクトルｃが、可変、すなわち時々刻々変化するものであっても良い。その場合、有効セット双対解記憶部４に保存する、有効セットＩに対する双対解候補ｙＩ、分配行列ＡＩ －１は、定数ではなく関数という形で保存することになる。具体的な例を以下に示す。制約係数行列Ａおよび費用ベクトルｃが、それぞれ次式で表現されているとする。なお、簡単のため、費用ベクトルｃは、変化せず、制約係数行列Ａが変化する場合を示す。
Figure 0006991326000009
 ・・・（９）
Figure 0006991326000010
 ・・・・・（１０）
ここで、制約次元数Mは、６、変数次元数Nは、８である。例えば、制御対象を宇宙機のスラスタとすると、ｕ，ｖ，ｗは、各スラスタにより生成される推力方向を表し、ｘ，ｙ，ｚは、各スラスタの位置を表し、式（９）の上３行は、各スラスタにより生成される推力方向を表し、下３行は、各スラスタの位置と各スラスタの推力方向との外積を表す。制約次元数Mが６、変数次元数Nが８であるから、双対解候補探索部で得られる情報として、すべての有効セットＩは、以下のようになる。
Figure 0006991326000011
 ・・・・（１１）
また、双対解候補探索部で得られる情報として、双対解候補ｙＩは、以下の形のように上記各スラスタにより生成される推力方向の変数のｕ，ｖ，ｗ，各スラスタの位置の変数ｘ，ｙ，ｚの関数として表される（一部抜粋）。
Figure 0006991326000012
 ・・・（１２）
式（９）～（１１）から、有効セットＩごとの双対解候補ｙＩは、式（１２）のように文字式によって表現される。解候補算出装置２は、文字式の形式で有効セット双対解記憶部４に、式（１２）のように文字式によって表現された双対解候補ｙＩを有効セットＩに関係づけて記憶する。
解候補算出装置２のＳｔｅｐ３は、制約係数行列Ａの第ｉ列の制約次元数MのM次元ベクトルをｉ列制約係数ベクトル ａｉ として、ｉ列制約係数ベクトルａｋ （ｋ∈Ｉ）をｋの小さい順に横に並べた制約次元数MのM×M行列を基底行列 ＡＩにする。この際、上記文字列で表現された制約係数行列Ａから、同じく文字列で表現された、制約次元数MのM×M行列を基底行列 ＡＩを求める。
さらに、基底行列 ＡＩ のランクが、制約次元数M、すなわちフルランクであるかどうかを確認する。なお、上記の文字列で表現された基底行列 ＡＩのランクは、文字列に依存した形で表現できるので、数式処理によってランクを求める。フルランクである場合は、Ｓｔｅｐ４へ進み、フルランクでない場合は、有効セットを記憶することなく、解候補算出装置２のＳｔｅｐ１へ戻り次の有効セットに移行する。
次に、最適解算出装置３の処理について説明する。最適解算出装置３のＳｔｅｐ１は、有効セット双対解記憶部４に記憶される、すべての有効セットＩに対して、制約ベクトルｂと双対解候補ｙＩとの内積を計算する。この際、最適解算出装置３は、有効セット双対解記憶部４に記憶される有効セットＩに関連付けられた双対解候補ｙＩを読出し、制約ベクトルｂと双対解候補ｙＩとの内積を計算する。
最適解算出装置３のＳｔｅｐ２は、Ｓｔｅｐ１で求めた内積が、最大となる有効セットＩを最適有効セットＩ＊ と定義し、式（７）によって、最適有効セットＩ＊を求める。
最適解算出装置３のＳｔｅｐ３は、標準形の線形計画問題の式（２）について、最適有効セットＩ＊を基底とした実行可能基底解を最適解ｘ＊ として求める。実行可能基底解の算出には、式（８）の線形連立方程式を解くことになる。この際、最適解算出装置３は、有効セット双対解記憶部４に記憶される最適有効セットＩ＊について、制約次元数Ｍ×Ｍの制約係数行列Aと制約次元数Ｍ次の列ベクトルである最適解ｘ＊ との積が、制約ベクトルｂとなることを表す主問題の制約条件式（２）を設定した連立方程式を解く。この際、制約係数行列Aの各係数に当たる関数（例えば、式（９））の各変数（式（９）の場合、ｕ，ｖ，ｗ，ｘ，ｙ，ｚ）に、実際の数値（例えば、上記の例では、各スラスタにより生成される推力方向、位置）を代入して制約次元数Ｍ×Ｍの制約係数行列Aを求めて、上記連立方程式を解き、最適解ｘ＊を得る。
なお、解候補算出装置２は、実施の形態１に記載の図１２と同様に、プロセッサ１００２、メモリ１００３、入出力インターフェース１００５のハードウェア構成で実現できる。最適解算出装置３は、実施の形態１に記載の図１３と同様に、プロセッサ１１０２、メモリ１１０３、入出力インターフェース１１０５のハードウェア構成で実現できる。
上述のように本実施の形態によれば、文字式関数として有効セット双対解記憶部４に保存する構成となっているため、最適解算出装置３で実行する前に、文字式に実際の値を代入して関数を計算し、制約係数行列Ａおよび費用ベクトルｃを求めて、解候補算出装置２の各Ｓｔｅｐを実行することで、制約係数行列Ａおよび費用ベクトルｃが固定でない場合にも実施の形態１，２と同様の手続きで主問題の最適解を反復なしで高速に算出することができる。
これは、例えば、仕様が異なるが同数のスラスタを有する宇宙機がある場合に、制約係数行列を文字式で表現することで、文字式の値に各仕様の値を入力するだけで、特定の仕様の宇宙機向けの解候補算出装置２を構成することができる。
実施の形態４．
本実施の形態は、上記実施の形態１ないし３の概念を用いて、複数のスラスタによる並進方向の位置および姿勢の制御を行う複数のスラスタを有する宇宙機についての実施の形態である。なお、複数のスラスタの数は、制御を行う自由度の数より多いとする。例えば、６自由度の制御を行う場合は、７本以上のスラスタを用いる。
図６は、複数のスラスタによる並進及び姿勢制御を行う宇宙機の構成を示すブロック図である。図において、宇宙機４００は、宇宙機４００を並進させる量および姿勢を算出する並進姿勢制御装置４０１、並進姿勢制御装置４０１からの指令を受けて、宇宙機４００のスラスタを制御するスラスタ制御装置４０２およびスラスタ制御装置４０２の制御対象となるスラスタ４０３を含む構成を有し、宇宙機の位置および姿勢を制御する。
スラスタ制御装置４０２は、並進姿勢制御装置４０１からの指令を受け、スラスタ４０３が出力すべき並進力およびトルクが入力される制御入力部４１１、制御入力部４１１からの信号（目標の並進力及びトルク）を受けてスラスタ４０３ごとに推薬を噴射する制御信号を演算するスラスタ駆動演算部４１２、およびスラスタ駆動演算部４１２の指示によって制御入力部４１１からの信号に基づきスラスタ４０３の噴射量最適化問題を上述の線形計画問題に変換した問題を解く最適解算出装置（回路）４１３を備える。
なお、制御入力部４１１は、数値的な安定性を向上させるため、力とトルクとでオーダーを合わせて規格化しても良い。また、制御入力部４１１は、無くても良く、並進姿勢制御装置４０１から入力された信号である、スラスタ４０３が出力すべき並進力およびトルクを直接スラスタ駆動演算部４１２に入力しても良い。
最適解算出装置（回路）４１３は、上述の実施の形態の有効セット双対解記憶部に対応する有効セット双対解記憶部４１４、スラスタ駆動演算部４１２から時々刻々変化する制約ベクトルｂに相当する制御ベクトルｂｓを受け取る制御ベクトル入力部４１５、有効セット双対解記憶部４１４から双対解候補ｙＩを読み出して、制御ベクトル入力部４１５から制御ベクトルｂｓを受け取り、制御ベクトルｂｓと双対解候補ｙＩとから最適有効セットＩ＊を求める最適有効セット演算部４１６、および最適有効セットＩ＊から最適解ｘ＊を求める最適解演算部４１７を備える。
以下、宇宙機４００の処理の流れについて説明する。
前準備として、ステップ４－０は、上記実施の形態の解候補算出装置２を用いて、有効セット双対解記憶部に有効セットＩと関連付けられた双対解候補ｙＩを記憶させる。
ここで、宇宙機４００のスラスタ４０３の総噴射量を最小化する問題は、標準形の線形計画問題である式（１）～（３）に帰着されることができるため、上記実施の形態を適用できる。以下、対応関係を説明する。
上記実施の形態の式（１）～（３）において、変数ベクトルｘは、最適化すべきスラスタ４０３の噴射量をスラスタ４０３の個数次元の非負の変数ベクトルとし、制約ベクトルｂは、上記目標並進力および目標トルクを縦に並べたベクトルとし、費用ベクトルｃは、最小化すべきもの（総噴射量）と変数ベクトル（各スラスタの噴射量）が定義されると決まる、すべての要素が１のベクトルとする。さらに、制約係数行列Ａは、スラスタ４０３の配置に依存する行列であり、式（９）の説明で示した、各スラスタによって生成される推力方向を表す行および各スラスタの位置と推力方向との外積を表す行で構成される。
ここで、変数ベクトルｘの次元数Ｎは、制約ベクトルｂの次元数Ｍよりも、大きいことを前提とし、上述の実施の形態１の式（１）から（３）で示した条件を満たす。これは、宇宙機を制御する自由度よりも、変数ベクトルｘの次元数Ｎ、すなわちスラスタ４０３の本数が多いことを前提としている。例えば、並進３自由度、回転３自由度を合わせて、６自由度の制御をする場合は、宇宙機が、７本以上のスラスタ４０３を有することを前提とする。宇宙機を制御する自由度数は、６に限らない。例えば、一部のスラスタ４０３が故障した場合、宇宙機を制御する自由度を５以下にして、制御することもできる。
以上のように構成することで、制御ベクトルｂ、費用ベクトルｃ、制約係数行列Aを定めた標準形の線形計画問題を解いて変数ベクトルｘを求めることは、上記目標推進力および目標トルクを満たしながら、総噴射量を最小にするスラスタ４０３の噴射量を求めることになる。
ステップ４－０は、上述の実施の形態の解候補算出装置２のいずれかを適用して、予め問題入力部２１が、スラスタ４０３の宇宙機４００における配置から制約係数行列Ａを求め、要素が１のＮ次元の費用ベクトルｃを入力し、双対解候補探索部２２が、上記宇宙機用の標準形の線形計画問題（式（１）から（３））の双対問題（式（４）、（５））を求め、求めた双対問題の有効セットＩと双対解候補 ｙＩを求めて、関連付けて、宇宙機４００の有効セット双対解記憶部４１４に記憶する。
具体的には、スラスタ４０３の宇宙機４００における配置から求めた制約係数行列Ａおよび要素が１のＮ次元の費用ベクトルｃで定義される標準形の線形計画問題の双対問題に対して、上記図２の処理フローと同じの処理フローを実行する。
また、実施の形態２のように制約係数行列Ａの逆行列である分配行列ＡＩ －１を算出し、有効セットＩごとに関連付けて、有効セット双対解記憶部４に分配行列ＡＩ －１を保存しても良い。また、実施の形態３のように各スラスタ４０３により生成される推力方向を表す上３行および各スラスタ４０３の位置と上記推力方向との外積を表す下３行をスラスタ４０３の推力方向、位置の関数によって制約係数行列Ａおよび費用ベクトルｃが変化するようにしても良い。
なお、費用ベクトルｃは、一部の成分だけ大きな値にすることで、対応するスラスタ４０３の噴射量が、ほぼ０となる。この特性を利用して、一部のスラスタ４０３が故障したときに、対応する費用ベクトルｃの成分を大きな値にして、近似的に同じ処理で計算することができる。
なお、ステップ４－０は、上述の実施の形態の解候補算出装置２で実行され、解候補算出装置２は、図１２の実施の形態１に記載の図１２と同様に、プロセッサ１００２、メモリ１００３、入出力インターフェース１００５のハードウェア構成で実現できる。
ステップ４－１は、並進姿勢制御装置４０１が、ＧＮＳＳ（Ｇｌｏｂａｌ Ｎａｖｉｇａｔｉｏｎ Ｓａｔｅｌｌｉｔｅ ｓｙｓｔｅｍ）例えば、ＧＰＳを用いて、並進位置を、スターセンサーを用いて、姿勢を計測することによって、宇宙機４００の自己位置、速度、姿勢、および角速度を求める。また、これらの目標値である宇宙機４００の目標位置、目標速度、目標姿勢、および目標角速度をもとに、並進及び姿勢制御に必要となる所望の並進力およびトルクである目標並進力および目標トルクを算出する。これは例えば、目標位置との差をフィードバックするＰＩＤ制御などの制御則により実現する。
ステップ４－２は、スラスタ駆動演算部４１２が、ステップ４－１で得られた目標並進力および目標トルクを実現するためのスラスタ４０３の噴射量を算出する。
上述の実施の形態の最適解算出装置３のいずれかを適用することによって、最適なスラスタ４０３の噴射量を高速に算出することができる。
具体的には、上記最適解算出装置３にあたる宇宙機のスラスタ制御装置の最適解算出装置４１３が、時々刻々変化する上記目標並進力および目標トルクを縦に並べたベクトルである制約ベクトルｂが、時間ごとに繰り返し与えられて、制約ベクトルｂに応じた最適解ｘ＊を高速に算出する。この処理フローは、図３と同様になる。
最適解算出装置４１３は、Ｓｔｅｐ１で、有効セット双対解記憶部４１４に記憶されるすべての有効セットＩに対して、上記目標並進力および目標トルクを含む制約ベクトルｂと双対解候補ｙＩとの内積を計算する。
最適解算出装置４１３は、Ｓｔｅｐ２で、Ｓｔｅｐ１で求めた内積が、最大となる有効セットＩを最適有効セットＩ＊ と定義し、式（７）によって、最適有効セットＩ＊を求める。
最適解算出装置４１３は、Ｓｔｅｐ３で、標準形の線形計画問題の式（２）について、最適有効セットＩ＊を基底とした実行可能基底解を最適解ｘ＊ として求める。なお、実行可能基底解の算出は、式（８）の線形連立方程式を解く。ここで、式（８）は、各スラスタ４０３によって生成される推力方向を表す行および各スラスタの位置と推力方向との外積を表す行で構成される制約次元数Ｍ×Ｍの制約係数行列Aと、最適化すべきスラスタ４０３の噴射量をスラスタ４０３の個数次元の非負の列ベクトルである最適解ｘ＊ との積が、上記目標並進力および目標トルクを含む制約ベクトルｂとなることを表す主問題の制約条件式（２）を設定した連立方程式である。
ステップ４－３は、スラスタ駆動演算部４１２が、ステップ４－２で最適解算出装置４１３が算出したスラスタ４０３の噴射量を含む指令をスラスタ４０３に送る。これによって、スラスタ４０３は、ステップ４－２で求めた噴射量を噴射し、総噴射量を最小にする宇宙機４００に所望の並進力およびトルクを与えることができる。
なお、スラスタ制御装置４０２は、実施の形態１に記載の図１３と同様に、プロセッサ１１０２、メモリ１１０３、入出力インターフェース１１０５のハードウェア構成で実現できる。入出力インターフェース１１０５は、制御入力部４１１の入力、スラスタ駆動演算部４１２からスラスタ４０３への出力を行い、プロセッサ１１０２は、スラスタ駆動演算部４１２、制御ベクトル入力部４１５、最適有効セット演算部４１６および最適解演算部４１７の演算を行い、メモリ１１０３は、スラスタ制御装置４０２の処理プログラムを保持し、有効セット双対解記憶部４１４を実現する。また、図１３の受信装置１１０６は、本実施の形態では、並進姿勢制御装置４０１からの情報を受信し、アクチュエータ１１０８は、スラスタ４０３にあたり、アクチュエータ制御装置１１０７は、スラスタ駆動演算部４１２の一部にあたる。
以上のように、本実施の形態によれば、並進および姿勢制御のアクチュエータとして複数のスラスタ４０３を有する宇宙機において、スラスタ４０３の配置から決定される行列Ａ、およびすべての要素が１である費用ベクトルｃを入力として、解候補算出装置２に動作させることにより得られる双対解候補および有効セットを有効セット双対解記憶部４に予め保存しておき、目標並進力および目標トルクを縦に並べた制御ベクトルｂを入力として最適解算出装置３を動作させてスラスタ４０３の噴射量を算出するから、総スラスタ噴射量を最小にするスラスタ４０３の噴射量を高速に算出することができる。
また、制御サイクルごとに制約ベクトルｂである目標並進力および目標トルクが変化するが、収束演算をすることなく、制御サイクル内の一定時間内に、制約を満たしながら総噴射量が最小となるスラスタ４０３の噴射量をリアルタイムで計算できる。宇宙機の搭載計算機は、信頼性や耐放射線性などが求められており、地上の計算機に比べて一般に低性能である。よって、最適な噴射量を一定時間で確実に算出できることは、信頼性を保証する上で重要である。本実施の形態によれば、宇宙機４００で変化しないスラスタ４０３の配置から決定される行列Ａ、およびすべての要素が１である費用ベクトルｃから求める双対解候補および有効セットを有効セット双対解記憶部４に予め保存するから、宇宙機４００のスラスタ制御装置で収束計算することなく、目標併進力および目標トルクを満たすスラスタ４０３の総噴射を最小にする噴射量を一定時間内に確実に算出できる。
実施の形態５．
実施の形態１では、非退化な標準形の線形計画問題を対象としていたが、本実施の形態では、標準形線形計画問題の双対問題を主問題とする非退化な線形計画問題を求解する線形計画問題求解システム、解候補算出装置、最適解算出装置について述べる。本実施の形態では、標準形線形計画問題の双対問題を主問題とする非退化な線形計画問題を対象とする。ここで、本実施の形態で対象とする標準形線形計画問題の双対問題を主問題とする非退化な線形計画問題は、与えられた制約係数行列Ａ及び制約ベクトルｃ並びに変数ベクトルｙの連立不等式で表される制約条件下で、変数ベクトルｙの関数である目的関数を最適化（最大に）する問題であり、次式で表される。
Figure 0006991326000013
 ・・・・・・（１３）
Figure 0006991326000014
 ・・・・・・（１４）
式（１４）において、ｃは、Ｎ次元の制約ベクトルであり、ｙは、制約ベクトルｃの次元数Ｎより小さいＭ次元の変数ベクトルであり、Ａは、変数ベクトルｙの次元数Ｍ×制約ベクトルｃの次元数Ｎの制約係数行列である。式（１３）において、ｂは、変数ベクトルｙの次元数Mの費用ベクトルであり、式（１４）を満たす変数ベクトルｙのうち、費用ベクトルｂの転置行列と変数ベクトルｙの内積を最大にするものを求めることを意味する。
式（１３）（１４）は、実施の形態１の式（１）～（３）で示される標準形の線形計画問題を主問題としたときの双対問題である式（４）（５）と同じ式であり、この点で、標準形線形計画問題の双対問題を主問題とする線形計画問題であるといえる。
なお、退化した問題である場合も摂動を加えることにより非退化な問題に変形することができるため、退化した問題も非退化な問題に変形して、上記問題として求解することができる。
本実施の形態は、式（１３）（１４）に示される線形計画問題の制約係数行列Ａと費用ベクトルｂは変化せず、制約ベクトルｃが変化する問題を繰り返し求解することを収束計算なしに高速に実現するものである。
図１は、本実施の形態の標準形線形計画問題の双対問題を主問題とする線形計画問題を求解する線形計画問題求解システムの構成を表す図でもある。実施の形態１と同様に、線形計画問題求解システム１は、入力に対して、線形制約を満たしつつ線形な評価関数を最小にする出力ベクトルを算出するものである。線形計画問題求解システム１は、解候補算出装置２及び最適解算出装置３により構成される。
解候補算出装置２は、対象とする線形計画問題の不変となる情報を入力として、入力された情報から対象とする線形計画問題を求め、求めた線形計画問題を主問題とする双対問題を求め、求めた双対問題の制約条件が有効となる実行可能な解候補と解候補に対応する有効セットをすべて求めて記憶部に記憶させる。
最適解算出装置３は、随時変化する変化情報を入力とし、解候補算出装置２で求めた解候補と入力の内積から最適となる有効セットを選択し、有効セットに対応する実行可能基底解を最適解として求めて出力する。
本実施の形態は、式（１３）（１４）に示される線形計画問題の制約係数行列Ａと費用ベクトルｂは変化せず、制約ベクトルｃが変化する問題を繰り返し求解することを収束計算なしに高速に実現するものである。以下、この処理の流れを説明する。
まず、解候補算出装置２は、有効セット双対解記憶部４に格納するための情報を算出する。これには、式（１３）（１４）で示される線形計画問題を主問題としたときの双対問題を求める。この双対問題は、式（１３）（１４）の式から下記のように表現できる。
Figure 0006991326000015
 ・・・・・・（１５）
Figure 0006991326000016
 ・・・・・・（１６）
Figure 0006991326000017
 ・・・・・・（１７）
式（１６）は、制約係数行列Ａと次元数Ｎの解ベクトルｘとの積が費用ベクトルｂとなる連立方程式を表す。式（１７）は、解ベクトルｘの要素が、０以上であることを表す。式（１５）は、式（１６）（１７）を満たす解ベクトルｘのうちで、制約ベクトルｃと解ベクトルｘの内積を最小にするものを求めることを意味する。
ここで、制約係数行列Ａと費用ベクトルｂが固定であるという前提から、式（１６）（１７）の制約条件を予め解くことができる。この双対問題の制約条件（式（１６）および（１７））を満足する領域（いわゆる実行可能領域）は、凸領域となる。一方、線形計画問題の一般的性質から、凸領域の頂点の一つが必ず最適解（式（１５）を最小とする解）となるといえる。すると、上記前提に基づけば、凸領域の頂点を双対問題の最適解候補として、予め全て求めておくことができる。
そこで、本実施の形態では、解候補算出装置２が、双対問題の最適解候補を全て求め、予め有効セット双対解記憶部４に記憶し、最適解算出装置３が、有効セット双対解記憶部４から最適解候補を読み出して式（１１）を最小にする最適解を求める構成をとる。
図７は、解候補算出装置２の処理のフローチャートを示す図である。以下、図に従って説明する。なお、解ベクトルｘが、次元数Ｎのベクトルであるため、双対問題の不等式制約（式（１６）（１７））は、Ｎ個の式となる。なお、以下の処理の前に、図１の双対解候補探索部２２の双対問題作成部２２１が、問題入力部２１から入力される標準形線形計画問題の双対問題を主問題とする線形計画問題から双対問題を作成する。
ステップ５－１は、解ベクトルｘが、次元数Ｎのベクトルであるため、主問題（標準形の線形計画問題の双対問題である式（１３）（１４））の制約ベクトルｃの次元数Ｎ個の不等式制約（式（１４））のうち、変数ベクトルｙの次元数Ｍ個の制約が有効となるとして、有効となる制約の組み合わせを有効セットＩと呼ぶ。このとき、相補性条件より、双対問題の不等式制約（式（１７））の内、残りの（Ｎ－Ｍ）個の制約が、有効になる。
ステップ５－１の処理としては、図１の双対解候補探索部２２の全有効セット作成部２２２が、変数ベクトルｙの次元数Ｍ個の制約が有効となる１つの組合せ（有効セットＩ）を設定する。なお、以下の他のステップから戻った場合は、新たな組合せを有効セットＩとして処理を続ける。
ステップ５－２は、全有効セット作成部２２２が、制約係数行列Ａから、有効セットＩに対応する基底行列となる行列ＡＩを作り、行列ＡＩのランクが、変数ベクトルｙの次元数M、すなわちするフルランクであるかを調べる。行列ＡＩがフルランクであれば、次のステップ５－３へ進み、フルランクでない場合は、ステップ５－１へ戻り、次の有効セットＩに移行する。
ここで、行列ＡＩは、制約係数行列Ａの第ｉ列の変数ベクトルｙの次元数Ｍのベクトルａｊ（１，，・・・・Ｎ）としたとき、ベクトルａｋ（ｋは、有効セットＩに属する）をｋの小さい順に横に並べた変数ベクトルｙの次元数M×Mの行列である。すなわち、行列ＡＩは、基底行列であり、このステップ５－２は、基底行列がフルランクであるかのチェックをすることになる。
ステップ５－３は、図１の双対解候補探索部２２の双対解算出部２２３が、次の連立方程式の解ｘＩ（M次元ベクトル）を算出する。
Figure 0006991326000018
 ・・・・・・（１８）
ステップ５－４は、双対解算出部２２３が、ステップ５－３で求めた解ｘＩの要素が全て正であるかを確認する。全て正である場合は、ステップ５－５に進み、要素の１つでも正でなければ、ステップ５－１に戻り、次の有効セットに移行する。なお、ステップ５－４において、制約係数行列Ａと制約ベクトルｃの値によっては、解ｘＩの要素に０が含まれる場合があるが、線形計画問題が非退化であるという前提から、解ｘＩの要素に０が含まれる場合は考慮する必要がない。
ステップ５－５は、双対解算出部２２３が、有効セットＩと対応するｘＩを双対解候補として、有効セットＩに関連付けて有効セット双対解記憶部４に記憶する。その後、全ての有効セットについて処理を終了していなければ、ステップ５－１に戻り、次の有効セットへ移行する。
次に、最適解算出装置３について述べる。図８は、最適解算出装置３の処理のフローチャートを示す図である。最適解算出装置３は、変化する制約ベクトルｃが繰り返し与えられて、それに応じた最適解ｙ＊を収束演算することなく高速に算出する。
ステップ６－１は、図１の最適解算出装置３の最適有効セット演算部３１が、有効セット双対解候補記憶部４に記憶された、すべての有効セットＩに対して、制御ベクトル入力部３０から入力された制約ベクトルｃＩと双対解候補ｘＩとの内積を計算する。ここで、ｃＩは、制約ベクトルｃの要素の中で、有効セットＩにおいて制約が有効としたｋ番目の制約式に対応する要素について、ｋの小さい順に並べたものであり、変数ベクトルｙの次元数ＭのＭ次元ベクトルである。
ステップ６－２は、最適有効セット演算部３１が、制約ベクトルｃＩと双対解候補ｘＩとの内積の中で内積が最小となる有効セットを求める。ここで、上記内積の中で内積が最小となる有効セットを最適有効セットI＊と呼ぶ。複数ある場合は、どれを選んでも主問題の最適解は変化しないため、その中の一つを任意に選ぶ。これを数式で表すと次式となる。
Figure 0006991326000019
 ・・・・・・（１９）
次に、ステップ６－３は、図１の最適解算出装置３の最適解演算部３２が、上記で求めた最適有効セットI＊について、主問題（式（１４））の制約条件を等式とした、制約係数行列Aの転置行列と解ベクトルｙとの積が制約ベクトルｃIとなる連立方程式を立て、この連立方手式の解を求める。この解が、求めるべき解ベクトルｙの最適解ｙ＊となる。この連立方程式は、次式で表される。
Figure 0006991326000020
 ・・・・・・（２０）
以上によって、制約ベクトルｃについて、主問題の解ベクトルｙの最適解ｙ＊を算出することができる。
なお、解候補算出装置２は、実施の形態１に記載の図１２と同様に、プロセッサ１００２、メモリ１００３、入出力インターフェース１００５のハードウェア構成で実現できる。また、最適解算出装置３は、実施の形態１に記載の図１３と同様に、プロセッサ１１０２、メモリ１１０３、入出力インターフェース１１０５のハードウェア構成で実現できる。
以上のように、標準形線形計画問題の双対問題を主問題とする線形計画問題を求解するのに、実施の形態１で説明した、求解対象の線形計画問題を主問題とする双対問題を求め、求めた双対問題の制約条件の有効セットを全て求めた上で、実行可能な双対解候補を求めて双対解候補として記憶し、双対解候補と時々刻々変化する制約ベクトルとの内積から最適有効セットを選択し、選択した最適有効セットに対応する実行可能基底解を求める手法を適用することになる。
この際、本実施の形態は、式（１３）（１４）の標準形線形計画問題の双対問題を求解対象の線形計画問題とする。この対象の線形計画問題は、式（１４）の随時変化する制約ベクトルｃの次元数である制約次元数Ｎは、変数ベクトルｙの変数次元数Ｍより大きい。また、制約条件の連立式は、式（１４）の連立不等式である。さらに実施の形態１の双対制約条件は、式（１６）の対象の線形計画問題の制約係数行列の転置行列で表される変数次元数Ｍ個の等式制約、および式（１７）の解ベクトルの非負制約の不等式制約である。また、最適解演算部の逆行列を求める行列は、制約係数行列Ａのうちの最適有効セットの識別子に対応する列ベクトルを除いた行列である。ここで最適有効セットの識別子とは、制約係数行列Ａのａｉのｉにあたるインデックスである。
本実施の形態によれば、標準線形計画問題の双対問題を主問題とする線形計画問題が与えられ、与えられた線形計画問題を主問題とする双対問題を求め、双対問題の制約条件が時間的に変化しないことを利用して、前記制約条件が有効となる組み合わせである有効セットＩごとに、前記双対問題の最適解候補をオフラインで算出し、前記有効セットとそれに対応する前記最適解候補を有効セット双対解記憶部に記憶したため、最適な有効セットである最適有効セットＩ＊、および対応する基底解を最適解ｙ＊として高速に算出することができる。
また、制約係数行列Ａ、費用ベクトルｂが固定であり、制約ベクトルｃが変化する非退化な標準線形計画問題の双対問題を主問題とする線形計画問題を対象とする線形計画問題求解システムにおいて、制約ベクトルｃに依存しない双対解候補ｘＩおよび前記双対解候補に対応する有効セットＩを探索し、有効セット双対解記憶部４に保存する解候補算出装置２と、有効セット双対解記憶部４に保存された双対解候補ｘＩと制約ベクトルｃの内積を計算し、最大となる有効セットＩ＊を基底とした主問題の実行可能基底解を最適解ｙ＊として算出する最適解算出装置３とを備えた構成となっている。このように構成したので、双対解候補ｘＩを予めオフラインで計算しておくことができ、時間ごとに変化する制約ベクトルｃが与えられてから、主問題の最適解ｙ＊を算出するまでの過程に反復計算が不要となる。このことによって、計算の高速化および信頼性の向上を実現することができる。また、最適解を算出するまでの時間を見積もることができるため、制御サイクル内に最適解を確実に得られ、最適性および制約条件を犠牲にすることなく実時間で制御することができる。
実施の形態６．
実施の形態６に係る線形計画問題求解システム１は、解候補算出装置２、最適解算出装置３、および有効セット双対解記憶部４から構成される点で、実施の形態５と同様であるが、実施の形態６では、図９に示すように、図７に示した解候補算出装置において有効セット双対解記憶部４に格納する情報を追加することで、さらに求解を高速化することができる。具体的には、解候補算出装置２のＳｔｅｐ５、および最適解算出装置３のＳｔｅｐ３を下記のように修正する。
解候補算出装置２のＳｔｅｐ５’は、有効セットＩに対応する双対基底解 ｘｂ Ｉを双対解候補ｘＩとして、有効セット双対解記憶部４に保存する（以上は、実施の形態１のＳｔｅｐ５に同じ）。さらに、有効セットＩに対する制約次元数Ｍ×Ｍの制約係数行列Ａの逆行列である分配行列ＡＩ －１の転置行列を算出し、有効セットＩごとに関連付けて、有効セット双対解記憶部４に分配行列ＡＩ －１の転置行列を保存する。有効セットＩに対する制約次元数Ｍ×Ｍの制約係数行列Ａの逆行列である分配行列ＡＩ －１を求めるには、例えば・・・とすれば良い。その後、解候補算出装置２のＳｔｅｐ１へ戻り、次の有効セットへ処理を移行する。
実施の形態５では双対解候補ｘＩを探索して終了したが、各最適解（有効セット）に対応した主問題の分配行列ＡＩ －１の転置行列もあらかじめ計算しておく。また、実施の形態５の最適解算出装置３において、図１０に示すように最適解算出装置３のＳｔｅｐ３を下記のように修正する。
最適解算出装置３のＳｔｅｐ３’は、有効セット双対解記憶部４に記憶されている最適有効セットＩ＊に関連付けられた、分配行列ＡＩ －１の転置行列を読出し、最適有効セットの分配行列ＡＩ* －１の転置行列とし、制約次元数Ｍ×Ｍの最適有効セットの分配行列とベクトルの積ＡＩ －１Tｂを求め、最適解ｙI* を計算する。すなわち、式（２０）の連立方程式を解く代わりに、行列とベクトルの積であるｙ＊＝ＡＩ －１Ｔを計算するから、計算量が少ない。Ｓｔｅｐ３’では、既に求まっている分配行列ＡＩ* －１の転置行列を用いるため、逆行列の計算をする必要がなく、最適解算出装置３の演算量をさらに少なくできる。
なお、解候補算出装置２は、実施の形態１に記載の図１２と同様に、プロセッサ１００２、メモリ１００３、入出力インターフェース１００５のハードウェア構成で実現できる。また、最適解算出装置３は、実施の形態１に記載の図１３と同様に、プロセッサ１１０２、メモリ１１０３、入出力インターフェース１１０５のハードウェア構成で実現できる。
本実施の形態の線形計画問題求解システム１によれば、解候補算出装置２で有効セットＩに対する制約次元数Ｍ×Ｍの制約係数行列Ａの逆行列である分配行列ＡＩ －１およびこの転置行列を算出し、有効セットＩごとに関連付けて有効セット双対解記憶部４に記憶し、最適解算出装置３が、有効セット双対解記憶部４に記憶されている最適有効セットＩ＊に関連付けられた、分配行列ＡＩ －１Tを読出して最適解ｙI* を計算するため、最適解算出装置３の演算量をさらに少なくできる。これによって、制御サイクルをさらに短くして、応答が早いまたは滑らかな制御を実現することができる。
実施の形態７．
本実施の形態は、実施の形態５または６の概念を用いて、複数のロータによる並進及び姿勢制御を行う複数のロータを有する飛翔体についての実施の形態である。この飛翔体は、ドローンであっても良い。推力を一方向にしか出力できないドローンは、推力の大きさよりも姿勢を変えるトルクを優先的に出力することで姿勢を迅速に変更させたい場合がある。ロータ回転数の上下限制約のもとで、可能な限り目標となる目標推進力と目標となる姿勢にする目標トルクを実現する問題は、目標並進力と目標トルクとに重み付き誤差を設けて、この重み付き誤差を最小にする線形計画問題に帰着できる。
なお、実施の形態４では、スラスタの数は、制御を行う自由度の数より多いとしたが、本実施の形態では、飛翔体の推力発生部の数に制限はない。
すると、上下限制約のあるロータ推力のもとで、目標となる並進力およびトルクとの重み付き誤差を最小化する問題は、線形計画問題を表す式（１３）（１４）を解き、解ベクトルｙを最適化することになるから、実施の形態５，６を適用することができる。以下、複数のロータに回転数の上下制約のあり、推力を一方向に出力する複数のロータを有するドローンを目標推力よりも目標姿勢を優先的に実現するロータ制御を行う装置について説明する。まず、実施の形態５の線形計画問題に、本実施の形態の問題を当て嵌める方法を説明する。
ドローン各ロータの推力を並べた列ベクトルをロータ推力ベクトルｕ, 各ロータが、ロータ推力ベクトルｕを発揮したときに、ドローンに出力される並進力およびトルクを並べた列ベクトルをドローン力ベクトルｆとする。すると、ドローン力ベクトルｆは、ドローンにおけるロータの配置及び各ロータの特性できまる定数行列Ｂを用いてロータ推力ベクトルｕとの関係が表される。（このあたりで列ベクトルの要素数を説明ください）すなわち、ドローン力ベクトルｆとロータ推力ベクトルｕとは、定数行列Ｂの各要素が定数であるから、下記に示す線形な関係となる。
Figure 0006991326000021
 ・・・・・・（２１）
ここで、定数行列Ｂは、ロータ配置及びロータ特性で決まる定数行列である。また、各ロータの推力であるロータ推力ベクトルｕには、上限制約ｕｍａｘおよび下限制約ｕｍｉｎがあり、これは以下の不等式で表される。
Figure 0006991326000022
 ・・・・・・（２２）
このとき、所望の目標となるドローンの目標並進力および目標トルクを並べた列ベクトルを目標ドローン力ベクトルｆｄ とすると、目標並進力及び目標トルクと、解となるロータ推力ベクトルｕによるドローン並進力及びトルクとの重み付き誤差を最小化する問題は、下記のように定式化できる。
Figure 0006991326000023
 ・・・・・・（２３）
Figure 0006991326000024
 ・・・・・（２４）
ここで、下添え字が１でベクトルが縦棒２本で挟んだ記号は、ベクトルの１ノルムを表す。式（２４）は、解となるロータ推力ベクトルｕの上下限制約を表す。式（２４）の条件のもと、式（２３）は、ドローン力ベクトルｆと目標ドローン力ベクトルｆｄとの差に、重み行列Wをかけたベクトルの１ノルムを最小にするものを求めることを表す。
上記式（２３）（２４）は、以下のようにすることで、標準形線形計画問題を主問題としたときの双対問題に対応する線形計画問題（実施の形態５の対象とする問題）と等価な線形計画問題に変形できる。
Figure 0006991326000025
 ・・・・・・（２５）
Figure 0006991326000026
 ・・・・・・（２６）
ここで、実施の形態５の制約係数行列の転置行列ＡＴ，費用ベクトルｂ，制約ベクトルｃは、下記のように表される。
Figure 0006991326000027
 ・・・・・・（２７）
Figure 0006991326000028
 ・・・・・・（２８）
Figure 0006991326000029
 ・・・・・・（２９）
ただし、ｎは、ドローンのロータ数、Ｉnは、ロータ数nのn次元単位ベクトル、 ０n×６は、ｎ×６の零行列、０ｎは、０を縦にロータ数ｎ個並べた零ベクトル、１６は、１を縦に６個並べたベクトルである。また、解ベクトルｙは、ロータ数ｎの要素を持つロータ推力ベクトルｕと、６要素を持つ仮想変数ηとを用いて y = [ｕＴ ηＴ]Ｔと表される。ここで、仮想変数ηの各要素は、ドローン力ベクトルｆと目標ドローン力ベクトルｆｄとの差に、重み行列WをかけたベクトルＷ（ｆ－ｆｄ）の各要素の絶対値を表し、また、数字の６は、並進位置の自由度３と回転自由度３の合計６自由度であるから６となる。
式（２６）において、ｃは、Ｎ次元の制約ベクトルであり、ｙは、制約ベクトルｃの次元数Ｎより小さいＭ次元の変数ベクトルであり、Ａは、変数ベクトルｙの次元数Ｍ×制約ベクトルｃの次元数Ｎの制約係数行列である。式（２５）において、ｂは、変数ベクトルｙの次元数Ｍの費用ベクトルであり、式（２５）は、式（２６）を満たす変数ベクトルｙのうち、費用ベクトルｂと変数ベクトルｙの内積を最大にするものを求めることを意味する。
式（２７）は、式（１４）の制約係数行列Ａの転置行列について、最初の n 行は、各ロータの推力が上限制約ｕｍａｘ以下であることを表す。次の n 行は、各ロータの推力が下限制約ｕmin以上であることを表す。残りの１２行は、式（２３）を線形計画問題に変換する際に加えられる不等式を表し、ドローン力ベクトルｆと目標ドローン力ベクトルｆｄとの差に、重み行列Wをかけたベクトルの各要素が、各要素の絶対値以下であり、かつ各要素の絶対値に－１をかけた値以上となる、当然に成り立つ関係式が仮想変数ηに関して表現されている。
式（２８）は、費用ベクトルｂについて、最初のロータ数ｎ個の要素が０であり、続く６個の要素が－１であることを表す。
式（２９）は、制約ベクトルｃが、最初のロータ数ｎ個の要素が、上限制約ｕｍａｘ、続くロータ数ｎ個の要素が、下限制約ｕｍｉｎの負ベクトル、続く仮想変数の数６個の要素が、重み行列Ｗと目標ドローン力ベクトルｆｄの積Ｗｆｄ、及び続く仮想変数の数６個の要素が、積Ｗｆｄの負ベクトル－Ｗｆｄであることを示す。
以上のようにすることで、ドローンの目標並進力及び目標トルクと、解となるロータ推力ベクトルｕによるドローン並進力及びトルクとの重み付き誤差を最小化する問題は、実施の形態５の式（１３）（１４）で表現される標準形線形計画問題の双対問題を主問題とする線形計画問題として表現できる。したがって、［ｕＴ ηＴ］Ｔで表現される解ベクトルｙの最適解ｙ＊を求めるのに、実施の形態５または６を適用することができる。最適解ｙ＊が求まれば、最適なロータ推力ベクトルｕ＊を取り出すことができる。以下、実施の形態５を適用した一例を示す。
なお、飛翔体の推力発生部、例えば、ドローンのロータの数に制限はないとしたが、これは、次の理由による。実施の形態４の場合には、のスラスタの本数Ｎが、変数ベクトルｘの次元数であり、変数ベクトルｘの次元数Ｎが、制約式の数Ｍより多いことが、前提であった。しかし、本実施の形態のドローンの場合には、制約式の数Ｎより変数ベクトルｙの次元数Ｍが多いことが前提となるところ、線形計画問題に変形した時点で、上記式（２８）（２９）において、ロータ数ｎの数に依らず、制約式の数Ｎより変数ベクトルｙの次元数Ｍが大きいことが常に成り立つからである。
具体的には、制約式の数Ｎは、ロータ推力の上限と下限の制約がロータ数ｎの２倍と、式（２３）を線形計画問題に変換する際に加えられる不等式の数１２とを加算した「２＊ｎ＋１２」であり、変数ベクトルｙの次元数Ｍは、ロータ数ｎと自由度数６とを加算した「ｎ＋６」であることから、必ず制約式の数Ｎが、変数ベクトルｙの次元数Ｍより大きくなるからである。
また、上記は、ドローンを６自由度制御する場合の例を示したが、通常のマルチコプターは、全てのロータ７０３が同じ方向を向いているから、並進１自由度と回転３自由度の４自由度しか制御できない。この場合でも、上記制約式の数Ｎは、「２＊ｎ＋８」、変数ベクトルｙの次元数Ｍは、「ｎ＋４」となるから、必ず制約式の数Ｎが、変数ベクトルｙの次元数Ｍより大きくなる。
図１は、本実施の形態の線形計画問題を求解する線形計画問題求解システムの構成を表す図でもある。実施の形態１と同様に、線形計画問題求解システム１は、入力に対して、線形制約を満たしつつ線形な評価関数を最小にする出力ベクトルを算出するものである。線形計画問題求解システム１は、解候補算出装置２及び最適解算出装置３により構成される。解候補算出装置２は、対象とする線形計画問題の不変となる情報を入力として、入力された情報から対象とする線形計画問題を求め、求めた線形計画問題を主問題とする双対問題を求め、求めた双対問題の制約条件が有効となる実行可能な解候補と解候補に対応する有効セットをすべて求めて記憶部に記憶させる。最適解算出装置３は、随時変化する変化情報を入力とし、解候補算出装置２で求めた解候補と入力の内積から最適となる有効セットを選択し、有効セットに対応する実行可能基底解を最適解として求めて出力する。
飛翔体制御装置は、変数次元数のロータを有する飛翔体制御装置であって、ロータ発揮すべき並進力およびトルクが入力される制御入力部と、並進力およびトルクの情報から制約ベクトルを作成して制約ベクトル入力部に入力し、最適有効セット演算部に最適有効セットを求めさせ、解演算部に最適有効セットから変数ベクトルの最適解を求めさせ、得られた変数ベクトルからロータに出力する信号を作成するロータ駆動演算部とを備え、制約次元数は、変数次元数より大きく、連立式は、連立不等式であり、双対制約条件は、制約係数行列の転置行列で表される変数次元数個の等式制約を含み、不等式は、解ベクトルの非負制約であり、最適有効セット演算部は、求めた内積が最小となる有効セットを最適有効セットとして求める。
図１１は、本実施の形態の複数のロータによる並進及び姿勢制御を行うドローンの構成を示すブロック図である。図において、ドローン７００は、ドローン７００を並進させる量および姿勢を算出する並進姿勢制御装置７０１、並進姿勢制御装置７０１からの指令を受けて、ドローン７００のロータを制御するロータ制御装置７０２およびロータ制御装置の制御対象となるロータ７０３を含む構成を有し、ドローンの位置および姿勢を制御する。
ロータ制御装置７０２は、並進姿勢制御装置７０１からの指令を受け、ロータ７０３が出力すべき並進力およびトルクが入力される制御入力部７１１、制御入力部７１１からの目標となる並進力およびトルクの信号を受けてロータ７０３ごとにロータ回転を制御する制御信号を演算するロータ駆動演算部７１２、およびロータ駆動演算部７１２の指示によって制御入力部７１１からの信号に基づきロータ７０３の回転量の最適化問題を上述の実施の形態５，６の線形計画問題に変換した問題を解く最適解算出装置（回路）７１３を備える。
なお、制御入力部７１１は、数値的な安定性を向上させるため、並進力とトルクとでオーダーを合わせて規格化しても良い。また、制御入力部７１１は、無くても良く、並進姿勢制御装置７０１から入力された信号である、ロータ７０３が出力すべき並進力およびトルクを直接ロータ駆動演算部７１２に入力しても良い。
最適解算出装置（回路）７１３は、上述の実施の形態５，６の有効セット双対解記憶部４に対応する有効セット双対解記憶部７１４、ロータ駆動演算部７１２から時々刻々変化する制約ベクトルｃに相当する目標ドローン力ベクトルｆｄを受け取る制御ベクトル入力部７１５、有効セット双対解記憶部７１４から双対解候補ｘＩを読み出して、制御ベクトル入力部７１５から制約ベクトルｃを受け取り、制約ベクトルｃと双対解候補ｘＩとから最適有効セットＩ＊を求める最適有効セット演算部７１６、および最適有効セットＩ＊から最適解ｙ＊を求める最適解演算部７１７を備える。
以下、ドローン７００の処理の流れについて説明する。本実施の形態では、ドローンの目標並進力と目標トルクとに重み付き誤差を設けて、この重み付き誤差を最小にする線形計画問題に帰着する。帰着される線形計画問題の式（２５）（２６）は、実施の形態５の式（１３）（１４）と同じであるから、結果的に実施の形態５を適用できる。
前準備として、ステップ６－０は、上記実施の形態５または６の解候補算出装置２を用いて、有効セット双対解記憶部に有効セットＩと関連付けられた双対解候補ｘＩを記憶させる。
ステップ６－０は、実施の形態５または６の解候補算出装置２を適用して、予め問題入力部２１が、ロータ７０３のドローン７００における配置等および重みから式（２７）により制約係数行列Aを求め、式（２８）から費用ベクトルｂを入力し、双対解候補探索部２２がドローン用の線形計画問題の双対問題を求め、求めた双対問題の有効セットＩと双対解候補 ｘＩを求めて、関連付けて、ドローン７００の有効セット双対解記憶部７１４に記憶する。
なお、ステップ６－０は、上述の実施の形態の解候補算出装置２で実行され、解候補算出装置２は、図１２の実施の形態１に記載の図１２と同様に、プロセッサ１００２、メモリ１００３、入出力インターフェース１００５のハードウェア構成で実現できる。
ステップ６－１は、並進姿勢制御装置７０１が、自己位置、速度、姿勢、および角速度、ならびにこれらの目標値である目標位置、目標速度、目標姿勢、および目標角速度をもとに、並進及び姿勢制御に必要となる所望並進力とトルクを並進姿勢制御部７０１において算出する。これは例えば、目標位置との差をフィードバックするＰＩＤ制御などの制御則により実現できる。
ステップ６－２は、ロータ駆動演算部７１２が、ステップ６－１で得られた目標並進力およびトルクを可能な限り実現するための各ロータ推力を算出する。このとき、上下限制約のあるロータ推力のもとで、目標となる並進力およびトルクとの重み付き誤差を最小化する問題は、上述のように式（２７）－（２９）により定義される線形計画問題式（２５）（２６）を解いて、各ロータの推力を含む最適解ｙ＊を求める。
この際、制御ベクトル入力部７１５が、ロータ駆動演算部７１２から目標ドローン力ベクトルｆｄを受け取り、式（２９）によって、目標ドローン力ベクトルｆｄ、重みベクトルW、上限制約ｕｍａｘおよび下限制約ｕｍｉｎから制約ベクトルｃを求める。最適有効セット演算部７１６が、有効セット双対解記憶部７１４に記憶された有効セットＩに対応する制約ベクトルｃIの転置行列と、双対解候補ｘIとの内積が最小の最適有効セットＩ、最適双対解ｘI*を求める。
ステップ６－３は、最適解演算部７１７が、ステップ６－２で求めた最適有効セットＩ＊について、主問題（式（２５）、（２６））の制約条件（式（２６））の一部を等式とした連立方程式（AＴ Ｉ＊ｙＩ＊ = cＩ＊。式（２０）と同様）を立て、この連立方手式の解を求める。この解が、求めるべき解ベクトルｙの最適解ｙ＊となる。この連立方程式は、式（２０）で表される。さらに、最適解ｙ＊は、［ｕＴ ηＴ］Ｔであるから、最適解ｙ＊から各ロータ推力ベクトルｕ＊を取り出す。
ステップ６－４は、ロータ７０３に、ステップ７－２で算出したロータ７０３の制御量を含む指令を送る。これによって、ロータ７０３は、ステップ７－２で求めたロータ７０３の制御量によって制御され、推力の大きさよりも姿勢を変えるトルクを優先的に出力することで姿勢を迅速に変更するドローン７００の制御をすることができる。
なお、ロータ制御装置７０２は、実施の形態１に記載の図１３と同様に、プロセッサ１１０２、メモリ１１０３、入出力インターフェース１１０５のハードウェア構成で実現できる。入出力インターフェース１１０５は、制御入力部７１１の入力、ロータ駆動演算部７１２からロータ７０３への出力を行い、プロセッサ１１０２は、ロータ駆動演算部７１２、制御ベクトル入力部７１５、最適有効セット演算部７１６および最適解演算部７１７の演算を行い、メモリ１１０３は、ロータ制御装置７０２の処理プログラムを保持し、有効セット双対解記憶部７１４を実現する。また、図１３の受信装置１１０６は、本実施の形態では並進姿勢制御装置７０１からの情報を受信し、アクチュエータ１１０８は、ロータ７０３にあたり、アクチュエータ制御装置１１０７は、ロータ駆動演算部７１２の一部にあたる。
以上のように、式（２５）（２６）の標準形線形計画問題の双対問題を主問題とする線形計画問題を求解するのに、実施の形態１で説明した、求解対象の線形計画問題を主問題とする双対問題を求め、求めた双対問題の制約条件の有効セットを全て求めた上で、実行可能な双対解候補を求めて双対解候補として記憶し、双対解候補と時々刻々変化する制約ベクトルとの内積から最適有効セットを選択し、選択した最適有効セットに対応する実行可能基底解を求める手法を適用することになる。
この際、本実施の形態は、式（２５）（２６）の標準形線形計画問題の双対問題を求解対象の線形計画問題とする。この対象の線形計画問題は、式（２９）の随時変化する制約ベクトルｃの次元数である制約次元数Ｎは、ロータ数ｎの要素を持つロータ推力ベクトルｕと、６要素を持つ仮想変数ηとを用いて y = [ｕＴ ηＴ]Ｔと表される変数（解）ベクトルｙの変数次元数Ｍより大きい。また、制約条件の連立式は、式（２６）の連立不等式である。さらに実施の形態１の双対制約条件は、式（１６）の対象の線形計画問題の制約係数行列の転置行列で表される変数次元数Ｍ個の等式制約、および式（１７）の解ベクトルの非負制約の不等式制約である。また、最適解演算部の逆行列を求める行列は、制約係数行列Ａのうちの最適有効セットの識別子に対応する列ベクトルを除いた行列である。ここで最適有効セットの識別子とは、制約係数行列Ａのａｉのｉにあたるインデックスである。
本実施の形態によれば、並進および姿勢制御のアクチュエータとして、ロータを持つドローンにおいて、ロータ配置，ロータ特性，重み行列から決定される行列ＡＴおよびすべての要素が０と－１で構成される費用ベクトルｂを入力として、解候補算出装置２によって得られる双対解候補ｙIおよび有効セットIを有効セット双対解記憶部４に保存し、各ロータ推力ベクトルｕの上下限制約、重み行列、所望の並進力及びトルクにより定まる制約ベクトルｃ を最適解算出装置３に入力することで、目標となる並進力およびトルクとの重み付き誤差を最小化するような上下限制約の範囲内に収まるロータ推力を高速に算出することができる。
なお、ロータの推力方向は、複数のロータ間で同じである必要はなく、様々な方向に推力を発生するようにロータが配置されていても良い。この場合、ロータの推力の向きは、行列ＡＴで指定する。ただし、ドローンの飛行中は、ロータの向きが変わらないものとする。
なお、ドローンのロータ配置によっては１方向にしか推力を出すことができない場合も多く、この場合は並進よりも姿勢の制御を優先するように重み行列を設定することで、目的にあった制御を行うことができる。本実施の形態によって、上下限制約を考慮した最適なロータ推力を高速に決定することができ、ドローンの墜落防止やアクロバティックな飛行を可能にする。
また、上記は、ドローンを中心に説明したが、複数のロータまたは推力発生器を有する飛翔体であっても良く、特にロータまたは推力発生器の推力方向が動的に変わらない飛翔体に本実施の形態を適用でき、同様の効果を有する。
１ 線形計画問題求解システム
２ 解候補算出装置
３ 最適解算出装置
４ 有効セット双対解記憶部
２１ 問題入力部
２２ 双対解候補探索部
２２１ 双対問題作成部
２２２ 全有効セット作成部
２２３ 双対解算出部
３１ 制御ベクトル入力部
３２ 最適有効セット演算部
３３ 最適解演算部
４０１ 並進姿勢制御装置
４０２ スラスタ制御装置
４０３ スラスタ
４１１ 制御入力部
４１２ スラスタ駆動演算部
４１３ 最適解算出装置
４１４ 有効セット双対解記憶部
４１５ 制御ベクトル入力部
４１６ 最適有効セット演算部
４１７ 最適解演算部
７０１ 並進姿勢制御装置
７０２ ロータ制御装置
７０３ ロータ
７１１ 制御入力部
７１２ ロータ駆動演算部
７１３ 最適解算出装置
７１４ 有効セット双対解記憶部
７１５ 制御ベクトル入力部
７１６ 最適有効セット演算部
７１７ 最適解演算部

Claims (9)

  1. 与えられた制約係数行列、随時変化する制約ベクトルおよび前記制約ベクトルの次元数である制約次元数と異なる変数次元数の変数を持つ変数ベクトルの連立式で表現される制約条件の下で前記変数ベクトルおよび費用ベクトルの関数である目的関数を最適化する線形計画問題の解を求める線形計画問題求解システムにおいて、
    前記制約係数行列および前記費用ベクトルの情報を入力として、前記線形計画問題を主問題とする双対問題および前記双対問題の制約条件の式のうち等号が成立する組合せを表す有効セットを全て求め、前記有効セットごとに前記制約条件が有効となる実行可能な双対解候補を求めて前記有効セットと関連付けて前記双対解候補を記憶部に記憶する双対解候補探索部と、
    前記制約ベクトルを入力とし、前記記憶部に記憶された前記双対解候補と前記制約ベクトルとの内積から最適となる前記有効セットである最適有効セットを選択し、選択した前記有効セットに対応する実行可能基底解を最適解として求めて出力する最適解算出装置とを備える線形計画問題求解システム。
  2. 前記双対解候補探索部は、前記双対問題が前記変数次元数個の式を持つ双対制約条件を含む問題として求め、
    前記有効セットは、前記双対制約条件の複数の前記式のうち前記変数次元数および前記制約次元数の少ない方の次元数個の等号が成立する組合せを前記式の識別子で表したものであり、
    前記双対解候補は、前記有効セットごとに前記識別子が表す前記双対制約条件の前記式を満たす解ベクトルであり、
    前記最適解算出装置は、前記制約係数行列のうち前記最適有効セットの情報を使って求めた行列の逆行列と前記制約ベクトルとの積から前記変数ベクトルの最適解を求める最適解演算部を備える請求項１に記載の線形計画問題求解システム。
  3. 前記制約次元数は、前記変数次元数より小さく、
    前記変数ベクトルは非負であり、
    前記連立式は、連立方程式であり、
    前記双対制約条件は、前記制約係数行列の転置行列で表される前記変数次元数個の不等式制約であり、
    前記最適解算出装置は、求めた前記内積が最大となる前記有効セットを前記最適有効セットとして求める最適有効セット演算部を含み、
    前記最適解演算部の逆行列を求める前記行列は、前記制約係数行列のうちの前記最適有効セットの前記識別子に対応する列ベクトルを並べた行列であることを特徴とする請求項２に記載の線形計画問題求解システム。
  4. 前記制約次元数は、前記変数次元数より大きく、
    前記連立式は、連立不等式であり、
    前記双対制約条件は、前記制約係数行列の転置行列で表される前記変数次元数個の等式制約および前記解ベクトルの非負制約の不等式制約を含み、
    前記最適解算出装置は、求めた前記内積が最小となる前記有効セットを前記最適有効セットとして求める最適有効セット演算部を含み、
    前記最適解演算部の逆行列を求める前記行列は、前記制約係数行列のうちの前記最適有効セットの前記識別子に対応する列ベクトルを除いた行列であることを特徴とする請求項２に記載の線形計画問題求解システム。
  5. 与えられた制約係数行列、随時変化する制約ベクトルおよび前記制約ベクトルの次元数である制約次元数と異なる変数次元数の変数を持つ変数ベクトルの連立式で表現される制約条件の下で前記変数ベクトルおよび費用ベクトルの関数である目的関数を最適化する線形計画問題の解を求める線形計画問題求解システムの解候補算出装置であって、
    前記制約係数行列および前記費用ベクトルの情報を入力として、前記線形計画問題を主問題とする双対問題および前記双対問題の制約条件である双対制約条件の式のうち等号が成立する組合せを表す有効セットを全て求め、前記有効セットごとに前記制約条件が有効となる実行可能な双対解候補を求めて出力する双対解候補探索部と、
    前記双対解候補探索部から出力される前記有効セットと関連付けて前記双対解候補を記憶部に記憶する記憶部とを備える線形計画問題求解システムの解候補算出装置。
  6. 与えられた制約係数行列、随時変化する制約ベクトルおよび前記制約ベクトルの次元数である制約次元数と異なる変数次元数の変数を持つ変数ベクトルの連立式で表現される制約条件の下で前記変数ベクトルおよび費用ベクトルの関数である目的関数を最適化する線形計画問題の解を求める線形計画問題求解システムの最適解算出装置であって、
    前記線形計画問題の情報を入力する問題入力部と、
    前記線形計画問題を主問題とする双対問題および前記双対問題の制約条件である双対制約条件の式のうち等号が成立する組合せを表す有効セットの情報と、前記有効セットごとに前記制約条件が有効となる実行可能な双対解候補とを関連付けて記憶する記憶部と、
    前記制約ベクトルの情報が入力される制約ベクトル入力部と、
    前記記憶部に記憶された前記双対解候補と前記制約ベクトルとの内積から最適となる前記有効セットを最適有効セットとして選択する最適有効セット演算部と、
    前記最適有効セット演算部で選択した前記有効セットに対応する実行可能基底解を最適解として求めて出力する最適解演算部とを備える線形計画問題求解システムの最適解算出装置。
  7. 請求項６に記載の線形計画問題求解システムの最適解算出装置を備え、前記変数次元数のスラスタを有する宇宙機のスラスタ制御装置であって、
    前記宇宙機の発揮すべき並進力およびトルクが入力される制御入力部と、
    前記並進力および前記トルクの情報から前記制約ベクトルを作成して前記制約ベクトル入力部に入力し、前記制約ベクトルを含む前記線形計画問題を前記問題入力部に入力し、前記最適有効セット演算部に前記最適有効セットを求めさせ、最適解演算部に前記最適有効セットから前記変数ベクトルの最適解を求めさせ、得られた前記変数ベクトルから前記スラスタに出力する信号を作成するスラスタ駆動演算部とを備え、
    前記制約次元数は、前記変数次元数より小さく、
    前記変数ベクトルは非負であり、
    前記連立式は、連立方程式であり、
    前記双対制約条件は、前記制約係数行列の転置行列で表される前記変数次元数個の不等式制約であり、
    前記最適有効セット演算部は、求めた前記内積が最大となる前記有効セットを前記最適有効セットとして求めることを特徴とする宇宙機のスラスタ制御装置。
  8. 請求項６に記載の線形計画問題求解システムの最適解算出装置を備え、前記変数次元数のロータを有する飛翔体制御装置であって、
    前記ロータが発揮すべき並進力およびトルクが入力される制御入力部と、
    前記並進力および前記トルクの情報から前記制約ベクトルを作成して前記制約ベクトル入力部に入力し、前記最適有効セット演算部に前記最適有効セットを求めさせ、最適解演算部に前記最適有効セットから前記変数ベクトルの最適解を求めさせ、得られた前記変数ベクトルから前記ロータに出力する信号を作成するロータ駆動演算部とを備え、
    前記制約次元数は、前記変数次元数より大きく、
    前記連立式は、連立不等式であり、
    前記双対制約条件は、前記制約係数行列の転置行列で表される前記変数次元数個の等式制約および前記双対問題の解ベクトルの非負制約であり、
    前記最適有効セット演算部は、求めた前記内積が最小となる前記有効セットを前記最適有効セットとして求めることを特徴とする飛翔体制御装置。
  9. 与えられた制約係数行列、随時変化する制約ベクトルおよび前記制約ベクトルの次元数である制約次元数と異なる変数次元数の変数を持つ変数ベクトルの連立式で表現される制約条件の下で前記変数ベクトルおよび費用ベクトルの関数である目的関数を最適化する線形計画問題の解を求める線形計画問題の求解方法において、
    前記線形計画問題の求解方法は、コンピュータにより実行され、
    前記制約係数行列および前記費用ベクトルの情報を入力として、前記線形計画問題を主問題とする双対問題および前記双対問題の制約条件の式のうち等号が成立する組合せを表す有効セットを全て求め、前記有効セットごとに前記制約条件が有効となる実行可能な双対解候補を求めて前記有効セットと関連付けて前記双対解候補を記憶部に記憶する双対解候補探索ステップと、
    前記制約ベクトルを入力とし、前記記憶部に記憶された前記双対解候補と前記制約ベクトルとの内積から最適となる前記有効セットを選択し、選択した前記有効セットに対応する実行可能基底解を最適解として求めて出力する最適解算出ステップとを備える線形計画問題の求解方法。
JP2020520941A 2018-05-23 2018-05-23 線形計画問題求解システム、解候補算出装置、最適解算出装置、宇宙機のスラスタ制御装置および飛翔体制御装置並びに線形計画問題の求解方法 Active JP6991326B2 (ja)

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
PCT/JP2018/019880 WO2019224954A1 (ja) 2018-05-23 2018-05-23 線形計画問題求解システム、解候補算出装置、最適解算出装置、宇宙機のスラスタ制御装置および飛翔体制御装置並びに線形計画問題の求解方法

Publications (2)

Publication Number Publication Date
JPWO2019224954A1 JPWO2019224954A1 (ja) 2021-01-14
JP6991326B2 true JP6991326B2 (ja) 2022-01-12

Family

ID=68615771

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2020520941A Active JP6991326B2 (ja) 2018-05-23 2018-05-23 線形計画問題求解システム、解候補算出装置、最適解算出装置、宇宙機のスラスタ制御装置および飛翔体制御装置並びに線形計画問題の求解方法

Country Status (4)

Country Link
US (1) US11958635B2 (ja)
EP (1) EP3798869A4 (ja)
JP (1) JP6991326B2 (ja)
WO (1) WO2019224954A1 (ja)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US11155356B2 (en) 2020-02-19 2021-10-26 Kitty Hawk Corporation Thrust allocation using optimization in a distributed flight control system
US11163273B2 (en) * 2020-03-02 2021-11-02 Mitsubishi Electric Research Laboratories, Inc. Active set based interior point optimization method for predictive control
CN111353121B (zh) * 2020-03-31 2023-04-11 中国空气动力研究与发展中心超高速空气动力研究所 一种用于确定航天器解体碎片不确定性参数分布的方法
JP7571428B2 (ja) 2020-09-09 2024-10-23 日本電気株式会社 情報処理装置、情報処理方法、および情報処理プログラム
WO2024214167A1 (ja) * 2023-04-11 2024-10-17 日本電気株式会社 最適化装置、最適化方法、および、最適化プログラム

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20120271528A1 (en) 2008-04-04 2012-10-25 Honeywell International Inc. Methods and systems for the design and implementation of optimal multivariable model predictive controllers for fast-sampling constrained dynamic systems

Family Cites Families (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP3390492B2 (ja) * 1993-07-16 2003-03-24 三菱電機株式会社 宇宙機の制御装置、およびそのスラスタ制御方法
DE19718922C1 (de) 1997-04-25 1998-04-02 Daimler Benz Ag Verfahren zur treibstoffminimalen, rechnergestützten Ansteuerung für beliebig an einem Raumfahrzeug angeordnete Düsen
JP2000142596A (ja) 1998-11-16 2000-05-23 Nec Aerospace Syst Ltd 人工衛星の姿勢制御方法及びその装置
JP2001184334A (ja) * 1999-12-24 2001-07-06 Nec Corp 線形計画問題の求解方法および線形計画問題求解装置
JP4156544B2 (ja) 2004-03-05 2008-09-24 三菱電機株式会社 出力分配算出装置
US8812421B1 (en) * 2008-06-09 2014-08-19 Euler Optimization, Inc. Method and apparatus for autonomous synchronous computing
US9418046B2 (en) * 2014-05-14 2016-08-16 International Business Machines Corporation Price-and-branch algorithm for mixed integer linear programming
JP6838353B2 (ja) * 2016-10-31 2021-03-03 日本製鉄株式会社 鋼材の山分け計画作成装置、鋼材の山分け計画作成方法、およびプログラム
US10740431B2 (en) * 2017-11-13 2020-08-11 Samsung Electronics Co., Ltd Apparatus and method of five dimensional (5D) video stabilization with camera and gyroscope fusion

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20120271528A1 (en) 2008-04-04 2012-10-25 Honeywell International Inc. Methods and systems for the design and implementation of optimal multivariable model predictive controllers for fast-sampling constrained dynamic systems

Also Published As

Publication number Publication date
EP3798869A4 (en) 2021-06-23
US11958635B2 (en) 2024-04-16
US20210094705A1 (en) 2021-04-01
EP3798869A1 (en) 2021-03-31
JPWO2019224954A1 (ja) 2021-01-14
WO2019224954A1 (ja) 2019-11-28

Similar Documents

Publication Publication Date Title
JP6991326B2 (ja) 線形計画問題求解システム、解候補算出装置、最適解算出装置、宇宙機のスラスタ制御装置および飛翔体制御装置並びに線形計画問題の求解方法
Ulrich et al. Nonlinear adaptive output feedback control of flexible-joint space manipulators with joint stiffness uncertainties
Vukobratovic et al. Applied dynamics and CAD of manipulation robots
Jiang et al. Improving low-thrust trajectory optimization by adjoint estimation with shape-based path
Gondelach et al. Hodographic-shaping method for low-thrust interplanetary trajectory design
CN112313672B (zh) 用于无模型强化学习的堆叠的卷积长短期记忆
JP7035734B2 (ja) 強化学習プログラム、強化学習方法、および強化学習装置
Zhu et al. Fault-tolerant control algorithm of the manned submarine with multi-thruster based on quantum-behaved particle swarm optimisation
Selfridge et al. A multivariable adaptive controller for a quadrotor with guaranteed matching conditions
Park Robust and optimal attitude control of spacecraft with disturbances
Ortega et al. New results on control by interconnection and energy-balancing passivity-based control of port-Hamiltonian systems
Liu et al. Optimal path planning of redundant free-floating revolute-jointed space manipulators with seven links
Sasaki et al. Convex optimization for power tracking of double-gimbal variable-speed control moment gyroscopes
Kim et al. Optimal actuator failure control using a homotopy method
Manh et al. An adaptive neural network-based controller for car driving simulators
Angeletti et al. Automated nested co-design framework for structural/control dynamics in flexible space systems
Kou et al. Constrained control allocation of a quadrotor-like autonomous underwater vehicle
Shi et al. Symbolic programming of a graph-theoretic approach to flexible multibody dynamics
Chamseddine et al. Trajectory planning/re-planning for satellite systems in rendezvous mission in the presence of actuator faults based on attainable efforts analysis
Jin et al. Fault tolerant attitude control for small satellites using single gimbal control moment gyros and magnetic torquers
Chen et al. A reconfiguration control scheme for a quadrotor helicopter via combined multiple models
Ding et al. Smooth and proximate time-optimal trajectory planning of robotic manipulators
Papakonstantinou et al. Global Steering for Control Moment Gyroscope Clusters Using Heuristic Variable Search Techniques
Vang et al. Geometric attitude control via contraction on manifolds with automatic gain selection
Kuseyri MIMO H∞ control of three-axis ship-mounted mobile antenna systems

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20200622

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20210706

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20210802

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20211109

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20211207

R150 Certificate of patent or registration of utility model

Ref document number: 6991326

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150