WO2007094451A1

WO2007094451A1 - 運動解析方法、運動解析装置、コンピュータプログラム、及び記録媒体

Info

Publication number: WO2007094451A1
Application number: PCT/JP2007/052845
Authority: WO
Inventors: Joshua Hale
Original assignee: Japan Science And Technology Agency
Priority date: 2006-02-16
Filing date: 2007-02-16
Publication date: 2007-08-23
Also published as: JPWO2007094451A1

Abstract

本発明の運動解析方法を実現する情報処理装置１００は、物体表面上の各点に印加された外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて前記係数行列から非特異小行列を抽出し、抽出した非特異小行列を用いて連立一次方程式の解を算出するＣＰＵ１０１を備える。定式化した連立一次方程式の解を算出する際、ある加速度間の拘束を無視することによって一次独立の接触問題として解を求め、解を求めた後に全ての加速度間の拘束をチェックする。

Description

明細書

運動解析方法、運動解析装置、コンピュータプログラム、及び記録媒体技術分野

[0001] 本発明は、行列方程式を用いて運動解析を行う運動解析方法、運動解析装置、コンピュータプログラム、及び記録媒体に関する。

背景技術

[0002] ヒューマノイドロボットの挙動を解析する動的シュミレーシヨンでは、外部から印加される外力（接触力および衝突)を注意深く取り扱う必要がある。近年、正確かつ有効な解を見出すために精力的な研究がなされており、ペナルティ手法、撃力ベースの手法、解析的手法などの手法により運動の解析が行われて、る。

[0003] ペナルティ手法は、拘束を違反した度合いに応じて力を加えることで、シミュレーシヨンを進めるうちに拘束違反を解消する手法である（例えば、非特許文献 1及び 2参照)。接触力を正確に求めたい場合などに用いられてきた。拘束力は、侵入量 (すなわち、拘束違反の量）と速度とから直接的に求まる力拘束が解消するまでシミュレ一シヨンを繰り返す必要があると、う問題点を有して、る。

[0004] 撃力ベースの手法では、瞬間的な撃力が物体の速度を変化させるとして衝突および接触を解析する手法である（例えば、非特許文献 3参照)。この手法では、物体同士の相互作用やエネルギーの伝搬に着目しており、シュミレーシヨンの正確さよりも、むしろ物体の挙動を視覚化することに主眼を置、て、る。

[0005] 一方、解析的手法は、幾つかの制約があるものの現実のロボットの挙動を正確にシミュレートするのに最も適している。解析的手法では、線形拘束条件を取り扱い、拘束条件の式と運動方程式とを連立させて解くことで、拘束条件を満たすために必要な力を計算する（例えば、特許文献 4及び 5参照)。このように解析的手法では、拘束条件を満たす力を解くため、シミュレーションの時間をある程度広くとっても正確なシミュレーシヨンを行うことができる。

非特許文献 1：ウィルヘルム (J. Wilhelms)、他 2名著、「動的アニメーション:相互作用と帘¹ J御 (Dynamic animation: interaction ana control)」、ビンュ 7ノレコンヒュ' ~~タ (Visu al Computer) , 1988年、第 4卷、 P283— 295

非特許文献 2 :パンネ（M. van de Panne)、他 1名著、「センサーァクチユエータネットワーク（Sensor- actuator network)」、プロシーディングスエーシーエムシーグラフ（

Proc. ACM SIGGRAPH) , 1993年、第 27 (2)卷、 Ρ335— 343

非特許文献 3 :シュミドル（H. Schmidl)、他 1名著、「剛体シュミレーシヨンのための瞬間的な接 fis (A fast impulsive contact suite for rigid body simulation)」、視覚ィ匕とコンピュータグラフィックス（Visualization and Computer Graphics)、 2004年、第 10卷、 P

189 - 197

非特許文献 4 :フェザーストーン（R. Featherstone)著、「ロボット力学のアルゴリズム（R obot Dynamics Algorithms)」、クノレワーァカァミック (Kluwer Academic)、 1987年非特許文献 5 :バラフ (D. Baraff)著、「非侵入剛体に対する接触力の計算 (Fast cont act force computation for nonpenetratingngid bodies;」、、ン' ~~グフノ 94 コンファレンス（SIGGRAPH 94 Conference) , 1994年、 Ρ23— 34

発明の開示

発明が解決しょうとする課題

[0006] し力しながら、前述した解析的手法では、全ての接触点について剛体の運動を表す行列方程式を解く必要があるため、接触点が多い場合などに多くの計算時間を必要とすることが問題となっていた。

[0007] また、剛体上の複数の点に作用する接触力と加速度との関係が一次独立となり、係数行列が正則でありさえすれば、 LU分解法、 QR分解法等の行列分解を用いることにより解を求めることが可能であるが、実際には、多くの接触状態において一次従属の拘束を含むことが知られており、係数行列について独立性を仮定することはできない。この場合、ある拘束条件を無視し、他のすべての拘束条件を満たす解を求める手法がよく採用されるが、接触力自体を求めるのと同様に複雑な計算を要するという問題点を有している。

[0008] 本発明は斯力る事情に鑑みてなされたものであり、接触点が多ぐかつ接触力とカロ速度との関係が一次従属となる場合であっても、各接触点に印加される外力及びカロ速度を速やかに求めることができる運動解析方法、運動解析装置、コンピュータプログラム、及び記録媒体を提供することを目的とする。

課題を解決するための手段

[0009] 第 1発明に係る運動解析方法は、表面上の複数の点に外力が印加された場合の物体の運動を解析する運動解析方法において、各点に印加された外力と、外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて前記係数行列力ゝら非特異小行列を抽出し、抽出した非特異小行列を用いて前記連立一次方程式の解を算出し、前記連立方程式の途中結果及び最終結果をメモリに記録することを特徴とする。

[0010] 本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしてヽるため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。

[0011] 第 2発明に係る運動解析方法は、前記行列分解は、 QR分解であることを特徴とする。

[0012] 本発明にあっては、行列分解として QR分解を用いるため、既知の演算によって非特異小行列が抽出される。

[0013] 第 3発明に係る運動解析方法は、 QR分解を用いて前記係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて 1又は複数回の列交換及びギブンス回転を行うことを特徴とする。

[0014] 本発明にあっては、 QR分解により係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルにつ、て列交換及びギブンス回転を行うようにしているため、計算コストが低減される。

[0015] 第 4発明に係る運動解析方法は、前記上三角行列から零である対角要素を探索し

、零である対角要素が探索された場合、該対角要素の右隣から順に非零である要素を探索し、非零である要素が探索された場合、該要素を含む列ベクトルと前記対角要素を含む列ベクトルとの間の列交換を行うことを特徴とする。

[0016] 本発明にあっては、上三角行列力零である対角要素と、その対角要素の右側の列ベクトルに含まれる非零の要素とを探索する。そして、両要素が探索された場合、探索された対角要素及び非零である要素を含む 2つの列ベクトルを交換するようにしているため、非特異小行列が速やかに抽出される。

[0017] 第 5発明に係る運動解析方法は、前記物体は、床面を歩行する歩行ロボットであり、前記外力は前記床面力受ける抗カであることを特徴とする。

[0018] 本発明にあっては、解析対象として歩行ロボットを採用し、歩行ロボットが床面から受ける抗カを接触力問題として解くことができる。

[0019] 第 6発明に係る運動解析装置は、表面上の複数の点に外力が印加された場合の物体の運動を解析する運動解析装置において、各点に印加された外力と外力が印カロされた場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化する手段と、行列分解を用いて前記係数行列力ゝら非特異小行列を抽出する手段と、抽出した非特異小行列を用いて前記連立一次方程式の解を算出する手段と、前記連立一次方程式の途中結果及び最終結果をメモリに記録する手段とを備えることを特徴とする。

[0020] 本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしてヽるため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。

[0021] 第 7発明に係るコンピュータプログラムは、コンピュータに、表面上の複数の点に外力が印加された場合の物体の運動を解析させるコンピュータプログラムにおいて、コンピュータに、各点に印加された外力と、外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化させるステップと、コンビユータに、行列分解を用いて前記係数行列力非特異小行列を抽出させるステップと、コンピュータに、抽出させた非特異小行列を用いて前記連立方程式の解を算出させるステップとを有することを特徴とする。

[0022] 本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしてヽるため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。 [0023] 第 8発明に係る記録媒体は、コンピュータに、表面上の複数の点に外力が印加された場合の物体の運動を解析させるコンピュータプログラムが記録されているコンビュータでの読み取りが可能な記録媒体において、コンピュータに、各点に印加された外力と、外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化させるステップと、コンピュータに、行列分解を用いて前記係数行列から非特異小行列を抽出させるステップと、コンピュータに、抽出させた非特異小行列を用いて前記連立方程式の解を算出させるステップとを有するコンビュータプログラムが記録されて、ることを特徴とする。

[0024] 本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしてヽるため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。

発明の効果

[0025] 第 1発明及び第 6発明による場合は、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしている。したがって、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除して線形独立の関係のみを残すことができ、物体に印加されている外力及び加速度を解析的に速やかに求めることができる。

[0026] 第 2発明による場合は、行列分解として特定の QR分解を用いるため、既知の演算によって非特異小行列を抽出することができる。

[0027] 第 3発明による場合は、 QR分解により係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルにつヽて列交換及びギブンス回転を行うようにして、るため、計算コストを低減することができる。

[0028] 第 4発明による場合は、上三角行列から零である対角要素と、その対角要素の右側の列ベクトルに含まれる非零の要素とを探索する。そして、両要素が探索された場合、探索された対角要素及び非零である要素を含む 2つの列ベクトルを交換するようにしているため、非特異小行列を速やかに抽出することができる。 [0029] 第 5発明による場合は、解析対象として歩行ロボットを採用し、歩行ロボットが床面力も受ける抗カを接触力問題として解くことができる。

[0030] 第 7発明及び第 8発明による場合は、上述した発明に記載の運動解析方法をコンピュータプログラム及び記録媒体に記録されたコンピュータプログラムにより実現できるため、適宜のコンピュータを用いて物体に印加されている外力及び加速度を解析的に速やかに求めることができる。

図面の簡単な説明

[0031] [図 1]本実施の形態に係る運動解析装置の内部構成を示すブロック図である。

[図 2]2つの手法の計算コストの比較結果を示すグラフである。

[図 3]列交換プロセスを用いた場合の計算コストを示すグラフである。

[図 4]消去及び挿入プロセスを用いた場合の計算コストを示すグラフである。

[図 5]非零要素の探索手法を模式的に説明する説明図である。

[図 6]本発明に係るコンピュータプログラムにより実行される処理の手順を説明するフローチャートである。符号の説明

100 情報処理装置

101 CPU

102 ROM

103 RAM

104 記憶装置

105 入出力 IF

106 キーボード

107 モニタ

108 補助記憶装置

109 バス

110 記録媒体

発明を実施するための最良の形態

図 1は本実施の形態に係る運動解析装置の内部構成を示すブロック図である。本実施の形態に係る運動解析装置は、パーソナルコンピュータ、ワークステーション等の情報処理装置 100により実現される。情報処理装置 100は、演算手段としての CP U101を備えており、この CPU101には、 ROM102、 RAM103、記憶装置 104、入出力 IF105、補助記憶装置 108等のハードウェアがバス 109を介して接続されている。

[0034] ROM102には、バス 109に接続された各種ハードウェアの動作を制御するための制御プログラムが格納されている。 CPU101は、この制御プログラムを RAM103上にロードして実行することにより、ハードウェア全体の動作を制御する。また、 RAMI 03には、運動解析の途中結果及び最終結果が適宜書き出される。

[0035] 記憶装置 104は、ハードディスクドライブを備えており、本実施の形態で説明する運動解析方法を具体的に実現するためのコンピュータプログラム、及びこのコンビユータプログラムを実行する際に必要となるデータが記憶される。

[0036] 入出力 IF105には、入力デバイスであるキーボード 106、及び出力デバイスであるモニタ 107が接続される。情報処理装置 100は、キーボード 106を通じて運動解析の開始指示、運動解析に必要なパラメータの入力を受付ける。また、情報処理装置 1 00は、キーボード 106を通じて入力されたパラメータ、前述したコンピュータプロダラムの演算結果等をモニタ 107上に表示する。

[0037] 補助記憶装置 108は、本発明に係るコンピュータプログラムが記録された FD、 CD

-ROM, DVD— ROM等の記録媒体 110からコンピュータプログラムを読み取るための FDドライブ、 CD— ROMドライブ、 DVD— ROMドライブ等を備えている。補助記憶装置 108により読み取られたコンピュータプログラムは記憶装置 104に記憶される。 CPU101は、必要に応じて記憶装置 104力も前述のコンピュータプログラムを R AM103上にロードして実行することにより、情報処理装置 100を本実施の形態に係る運動解析装置として動作させる。

[0038] 情報処理装置 100による運動解析の対象として歩行ロボットを採用し、歩行ロボットが床面を歩行する際の運動を解析する手法を以下で説明する。歩行ロボットは、床面からの抗カを受けるため、この抗カ（接触力）と加速度との関係を運動方程式により表す。典型的な接触力の問題は、次の式のように表すことができる。 [0039] [数 1] a = Af + a₀

[0040] ここで、 fは物体表面の 1又は複数の接触点における接触力を表すベクトルであり、 aは各接触点での法線方向及び接戦方向の加速度を表すベクトルである。 Aは接触力と加速度とを関連付ける行列であり、物体の質量に関係する量を表す。また、 a

0は物体に接触力が印加されて、な、場合の前記接触点での加速度を表すベクトルである。

[0041] 行列 Aおよびベクトル a の各要素については、例えば、 SDFASTのような動力学

0

のソフトウェアを用いることにより、解析的に又は実験的に導出することができる。すなわち、行列 Aおよびベクトル a の各要素の値は既知の量と考えることができるため、

0

数式 1は、行列 Aを係数行列とした加速度 aおよび接触力 fに関する連立一次方程式と見ることができ、この連立一次方程式を解く問題に帰着することができる。

[0042] 係数行列を用いて連立一次方程式を解く方法としては、ガウスの消去法、 LU分解法、 QR分解法等が知られている。このうち、 QR分解に基づく方法は、係数行列が正則でありさえすれば解を求めることが可能である。しかし、行列が特異に近い場合には、解の精度を高めようとすると軸列の選択が必要なことが知られており、軸選択を行わない並列解法には適用範囲に限界がある。

[0043] 実際、数式 1に示した運動方程式では、多くの接触状態において一次従属の拘束を含むことが知られており、各加速度について独立性を仮定することができない。この場合、ある拘束条件を無視し、全ての拘束条件に適合する解を求める手法がよく採用されているが、接触力自体を求めるのと同様に複雑な計算を要する。また、摩擦力に起因する加速度によっては、任意の加速度の任意の次元にぉ、て一次従属となり得る。

[0044] そこで、本実施の形態では、行列 Aから非特異小行列を構築する手法 (以下、干渉性利用法という）について説明する。この干渉性利用法では、行列 Aについて消去しないサブセットの任意のインデックスを受付ける。このアルゴリズムでは、ある加速度間の拘束を無視することによって一次独立の接触問題として解を求め、解を求めた後に全ての加速度間の拘束をチェックする。不要な拘束を全て不要セットに置き換えた場合、これらの不要セットのインデックスの 1つは保持セットに含まれ、行列 Aの異なる一次独立の小行列を算出することができる。どの拘束が必要である力、どの拘束が他の拘束に影響を及ぼすかと、うことにつ、ては一般的には分力な、ため、保持する拘束のどの組み合わせについても、全ての拘束を満たす解が見つ力るまで線形独立の小行列を計算する。あるいは、元の行列 Aが悪状態又は不十分なランクであるために完全な解を求めることが出来ない場合もあり、そのときは最も拘束を満足する解が見つ力るまで全ての小行列について計算を行う。

[0045] 特異小行列の抽出プロセスには QR分解の理解が不可欠であるため、まず、列交換を保持するプロセスについて説明する。 QR分解では、 m行 n列の任意の行列 Aを、 m行 m列の直交行列 Qと m行 n列の上三角行列 Rとに行列分解する。

[0046] [数 2]

A = QR

[0047] 任意の行列 Aを QR分解することは可能である力直交行列 Q及び上三角行列 Rは一意ではない。 QR分解は、元の行列 Aが調整される場合、又は最初から計算をし直す努力を避けたい場合において特に有用である。行または列が追加された場合、 Q R分解はアップデートされ、行または列が削除された場合、 QR分解をダウンデートされる。または、ランク 1の行列の加算又は減算に続いて再計算が実行される。

[0048] QR分解のアップデート及びダウンデートは、行列の任意の列の交換に続いて分解を保持するという標準的な手法が採用される。例えば、第 n列と第 n+k列との間の交換を行う場合、第 n列および第 n+k列を消去し、その新たな 2箇所に列を挿入し、それぞれ 2回のアップデート及びダウンデートを行う。しかしながら、これらの演算は単一列の交換演算により統合することができ、演算資源をかなり削減することができる。

[0049] QR分解のアップデート及びダウンデートはギブンス回転を用いて実行される。ギブンス行列 Gn(i, j, θ )は、 η次元空間のうちの 2つの部分空間（1つは i次元の部分空間、もう 1つは j次元の部分空間とする）の間の回転を表している。ギブンス行列 Gn(i , j, θ )は、 η Χ ηの単位行列のうち、 (i, i)、 (i, j)、 (j , i)、 (j, j)の各要素をそれぞれ cos 0、—sin 0、 sin 0、 cos Θにより置き換えた行列である。 QR分解は以下のように与えられる。このとき、 Q 'の直交性は保持されるが、 Rは上三角行列でなくなる。

[数 3]

A = QR QG^T GR

R'

Q'

[0051] 行列 Aの QR分解が与えられた場合、以下のようにして第 k列を消去して行列分解を保持することができる。

[0052] [数 4]

A = QR = Q _r T _r ^Τ r ^T r ^Τ r ^T

r0 ^rk-\ ^rk ^rk+\ ^rn

T T T T

r0 " '^rk-\ ^rk+\ " ' ^rn

[0053] 例えば、 6行 6列の行列 R力第 3列の削除を考えた場合、行列 Rの構造は以下のように影響を受ける。

[0054] [数 5]

ここで、「 X」のシンボル、及び「 # 1」〜「 # 3」のシンボルで表される各要素は一般に非零の値をとる（但し、偶然的に零の値をとり得る）。また、シンボルが付されていない各要素は零の値をとる。行列 R⁽¹⁾は現時点で上三角行列ではないが、ギブンス回転を 3回適用することによって対角成分の直下の要素（すなわち、「# 1」〜「# 3」のシンボルで表した各要素）を 1つずつ零にすることができ、行列 "の上三角性を修復することができる。ここで、ギブンス回転は、インデックスを付した順に隣り合う行同士について適用され、その結果、対角要素の下側に非零の要素が生成されないことが保証される。

[0056] 行列 Aに列を挿入し、 QR分解をアップデートする場合も同様である。行列 Aの QR 分解が与えられた場合、以下のようにして、行列 Aに新しい列 c^Tを追加し、行列分解を保持することができる。

[0057] [数 6]

A τ Ί

"0 …" 1 c a_k ■a

= eR' = e| '… —】^Γ Q^TC^T r -r„

[0058] 例えば、 5行 4列の新たに第 3列を挿入することを考えた場合、行列 Rの構造は以下のように影響を受ける。

[0059] [数 7]

[0060] 上の式にぉ、ても「 X」のシンボル、及び「 # 1」〜「 # 3」のシンボルで表される各要素は一般に非零の値をとる（但し、偶然的に零の値をとり得る）。また、「のシンポルは、元々は零の要素である力ギブンス回転を適用することによって非零の値に変化し得る要素を示している。例えば、「# 1」の要素を零にするためにギブンス回転を適用した場合、「· 1」の要素が非零の値になり、「# 2」の要素を零にするためにギブンス回転を適用した場合、「· 2」の要素が非零の値になる。したがって、この場合には行列 R⁽¹⁾を上三角行列に変換するために 2回のギブンス回転を行う必要がある。

[0061] 同様の手法は、行列 Rから 2つの列を入れ替えた結果として生成される行列 R⁽¹⁾の構造を修復する際にも適用される。 2つの列の入れ替えは、列の挿入と削除とを同時的に実行することと等価であり、 R⁽¹⁾の構造を修復する前にアップデートと組み合わせて複数回の列の消去と挿入とを実行する場合と比べた場合、演算資源を削減する効果が高い。行列 Aの第冽と第 j列との入れ替えを考えた場合、以下のように行列分解を保持することができる。

[0062] [数 8]

[0063] 行列 R⁽¹⁾の構造はもはや上三角行列ではない。 8行 8列の上三角行列を例にとり、その第 3列と第 6列とを入れ替えた行列 R⁽¹⁾の構造を修復するプロセスにつ、て以下に説明する。行列 R⁽¹⁾は以下のように表すことができる。

[0064] [数 9]

[0065] この例における行列 R⁽¹⁾は、上三角行列の第 3列と第 6列とを入れ替えた行列であるから、その第 3列には対角要素より下側に 3つの非零の要素が存在する。これらの要素を零にするためにギブンス回転を行う。このとき、ギブンス回転の影響を受ける要素の数を最小にするために、隣り合う行同士にギブンス回転を適用する。「# 1」のシンボルで示した要素を零にするためにギブンス回転を行った結果、「· 1」のシンポルで示した要素は零から非零の値に変わる。「· 2」、「· 3」のシンボルで示した要素についても同様である。

[0066] 変形した結果の行列を R⁽²⁾とした場合、行列 ³は以下のように表すことができる。

[0067] [数 10]

[0068] 上の式で示したように、行列 ³には、第 4列及び第 5列の対角要素の下側に「 # 4 」及び「# 5」のシンボルで示した非零の要素が存在する。この第 2段階での変形では、第 1段階の 3回のギブンス回転により生成されたこれらの非零の要素を零にするために、 2回のギブンス回転を行う。その結果、「# 4」及び「# 5」で示した要素は零となり、「·4」及び「· 5」で示した要素はそのときのギブンス回転の影響を受けて値が変化する。この第 2段階での変形を行った結果、前述の行列 R⁽²⁾は上三角行列 R⁽³⁾ に変形される。対応する直交行列 Q(^¾は、直交行列 Qに逆の回転演算を施すことにより得られる。

[0069] 次に、計算コストについて説明する。

M行 N列の行列 Rから第 n列を消去した場合、その構造を修復するために必要なギブンス回転の数は、 max{min{M—n, N—n}， 0}となる。次いで、 n番目の位置に列を挿入した場合、修復する際に必要なギブンス回転の数は、 max{M-n, 0}となる。したがって、 N行 N列の正方行列について、列の消去及び挿入を行うことにより第 n列と第 n+k列との入れ替えを行う場合、修復に要するギブンス回転の回数は、 4N — 2 (2n+k)—lとなる。

[0070] 一方、列交換の手法を用いる場合、前述したように、第 1段階で k回のギブンス回転を行い、第 2段階で k—l回のギブンス回転を行うため、修復するために要するギブンス回転の数は 2k—lでよい。すなわち、トータルで 4 (N— (n+k) )回のギブンス回転を節約することができる。交換すべき 2つの列のインデックスのうち、大きい方のインデックスにより節約できるギブンス回転の数が定まり、インデックスが元の行列の最後の列に近づくにつれて減少する。

[0071] し力しながら、ギブンス回転は非零の要素を含んだ列に対してのみ適用すればよいので、前述した計算は必ずしも真の節約量を反映していない。この点を考慮し、行列回転を用いるよりも列回転を用いて節約量を計算することができる。 N行 N列の行列について第 n列と第 n+k列との間の列交換を行う場合、以下の回数だけギブンス列回転が必要となる。

[0072] [数 11]

[0073] 同様に、 M行 N列の行列力列の消去を行う場合、以下の回数だけギブンス列回転が必要となる。但し、 Mは、 ME [N, N+ 1]を満足する整数である。

[0074] [数 12]

-(N - n)(N + l - n)

[0075] 更に、列の挿入を行う場合には、以下の回数だけギブンス列回転が必要となる。伹し、 Mは、 ME [N+ 1, N + 2]を満足する整数である。

[0076] [数 13] (M - N) + -(N - n + l)(N - n + 2)

[0077] したがって、 2つの列を消去して 2つの列を挿入する場合の計算コストは、以下のように表すことができる。

[0078] [数 14]

{N - n - kf + 2 (N - n - k) + (N - nf

[0079] 図 2は 2つの手法の計算コストの比較結果を示すグラフである。このグラフでは 24行 24列の行列について計算コストを算出しており、底面の 2つの軸はそれぞれ交換対象の列のインデックスを示し、高さ方向の軸は相対コストを示している。このグラフから、消去及び挿入プロセスに比べて列交換プロセスは大部分の条件下で計算コストが良いことが分かる。相対コストが 1となる条件をグラフ中の破線により示している。すなわち、交換対象の列のインデックスの何れか一方が 20以下である場合、計算コストは交換プロセスを使用する方が良ぐ交換対象の列のインデックスが元の行列の最後の列（第 24列）に近づくにつれて計算コストは消去及び挿入プロセスを使用する方が良くなる。

[0080] 図 3は列交換プロセスを用いた場合の計算コストを示すグラフであり、図 4は消去及び挿入プロセスを用いた場合の計算コストを示すグラフである。何れのグラフにつヽても 24行 24列の行列について計算コストを算出している。何れのプロセスにおいても計算コストは元の行列のサイズに依存する。行列サイズが大きくなるにつれ、列交換プロセスの計算コストは、消去及び挿入プロセスの計算コストの 50%に近づくことが分力つた。このファクタ一は、任意の 2つの列の交換の可能性を等確率と仮定し、列回転を実行する回数をカウントすることにより見積もった値である。 4行 4列の行列、又はそれ以上のサイズを有する行列では、列交換プロセスは消去及び挿入プロセスよりも計算速度が平均的に大きぐ 3行 3列の行列で計算コストが等しくなる。

[0081] 次に、最大の非特異小行列を抽出する手法について説明する。与えられた行列から最大限に大きな非特異小行列を抽出する手法は、前述したような QR分解を基にしており、以下の仮定を必要とする。まず、与えられた行列が対象行列、又は行交換した対象行列であることを仮定する。また、正規直交の QR分解において行列 Rが与えられていることを仮定する。更に、消去することができない行および列の組を示す保持用インデックスのセットが与えられてヽることを仮定する。更に、保持用インデックスのセットの軽微な変更 (典型的には、新たなインデックスの追加）に続いて小行列の抽出を行うことを仮定する。更に、一連の類似の行列について一力小行列の抽出を繰り返すことを仮定する。

[0082] 行列の左上側のブロックの非特異小行列と明確に区別するために元の行列の行および列を再配列することを行う。再配列するインデックスは、続くプロセスの繰り返しの間だけ内部的に保持される。

[0083] 前述の仮定を用いて非特異小行列を抽出するために、対象行列 Aに対する QR分解を以下のように表す。

[0084] [数 15]

[0085] ここで、 0又は 3のインデックスを付した小行列（例えば、 A 、 Α 等）は正方行列で

0 3

あり、同じインデックスを付した全ての小行列は同じ次元を有する。小行列 R の対角

0 要素の全てを非零と仮定したとき、小行列 A が非特異であり、かつ、小行列 A の外

0 0 部に行または列の組を追加することによって拡張した小行列 A は特異行列となるこ

0

とが証明できる。

[0086] QR分解を構築することは、行列 Aの各列に亘つて互いに直交するベクトルのセットを見つけることと等価である。また、 QR分解を区分ィ匕することは、基を構成するまで直交するベクトルのセットを縮小することと等価である。すなわち、行列 Rの左上隅から完全なシリーズとなるように、行列 Rの対角要素にできるだけ多くの非零の要素を配列させるプロセスが要求される。 [0087] そのため、行列 Rの対角要素を左上隅力順に調べ、零の要素を探索する。本実施の形態では、微小な値を閾値として予め設定しておき、その閾値より小さい値を見つけた場合、その要素は零であると判定する。零の要素を見つけた場合、その要素のすぐ右隣から列のインデックスが増える方向に非零の要素を探索する。その行に非零の要素が見つからない場合には、行を一段だけ下げて対角要素の右側を探索し、非零の要素が見つ力るまで行を一段ずつ下げながら行列 Rの右下隅に至るまで探索を繰り返す。非零要素が見つカゝらない場合、行列 Rは区分化された構造を有しているため、本プロセスを終了する。なお、対角要素の下側は零であることが保証されているため、探索を省略できることはいうまでもない。図 5は非零要素の探索手法を模式的に説明する説明図である。

[0088] このような要素が見つかった場合、その要素を含む列は零対角要素を含んだ列と交換され、アップデートされた行列 Rの次の対角要素が零であるか否かの判定が行われる。列交換は、非零対角要素の上側及び左側には影響を及ぼさないため、行列 Rの構造は保持される。小行列における零対角要素の下側及び右側に非零要素が存在するということは、その零対角要素を含んだ列と行列 Rの非零対角要素を含んだ各列とがー次従属でな!、ことを保証する。

[0089] 保持用セットは、行列 Rの構造を修復するのに先立ち、保持した列を非特異ブロックの最初の列と交換することにより適用することができる。例えば、ある列が保持用セットに存在するが非特異ブロックに存在しない場合、その列はそれ自身が保持用セットに含まれない非特異ブロックの最初の列と交換される。非特異ブロック中の目的の列の対角要素が零になるという結果である場合には、保持セットは一次従属の関係にある列を含み、もはや解が存在しないためプロセスを終了する。

[0090] この列交換の手続きは、保持用セットの操作は最小の列交換アップデートにより処理することができ、アルゴリズムの次の反復の間で期待された干渉性を使用することができることを確実にする。また、保持用セットが実現可能である場合には速やかに決定することができる。このパターンを用いてオリジナルのマトリックスの順序を入れ替えることにより、最初の QR分解は区分ィ匕した形態になる力、又は近づき、必要な列交換の数を減らすことができる。 [0091] このようにして求めた解は対象行列と同等に扱うことができる力結果として導出された非正則ブロックはアルゴリズムの実行を繰り返す度にチェックする必要がある。これは、候補となるブロック行列を完全に QR分解することを要求する。 Q行列がァルゴリズムを通して保持される場合、分解を再計算する必要はなくなるが、候補となる行列のオリジナルの行列に対する相対サイズに依存することとなり、 Q行列の間接的な保持及び QR分解を一力も再計算することを避けることが十分可能になる。

[0092] 次に、本発明に係るコンピュータプログラムにより実行される処理の手順を説明する。図 6は本発明に係るコンピュータプログラムにより実行される処理の手順を説明するフローチャートである。外部力も入力されたパラメータに基づいて物体 (例えば、歩行ロボット）の運動解析を行う場合、物体表面上の各点に印加された外力と加速度との関係を定式化する (ステップ Sl)。物体表面上の各点に印加された外力と加速度との関係は、前述の数式 1に示したように連立一次方程式を用いることにより、物体表面の 1又は複数の接触点における接触力を表すベクトル f、各接触点での法線方向及び接線方向の加速度を表すベクトル a、接触力と加速度とを関連付ける行列 Aにより記述することができる。

[0093] 次、で、行列分解 (特に、 QR分解)を用いることにより行列 Aから非特異小行列を抽出し (ステップ S 2)、抽出した非特異小行列を用いて定式化した連立一次方程式の解を算出する (ステップ S3)。このとき、ある加速度間の拘束を無視することによつて一次独立の接触問題として解を求め、解を求めた後に全ての加速度間の拘束をチェックする。不要な拘束を全て不要セットに置き換えた場合、これらの不要セットのインデックスの 1つは保持セットに含まれ、行列 Aの異なる一次独立の小行列を算出することができる。どの拘束が必要である力、どの拘束が他の拘束に影響を及ぼすかということについては一般的には分力もないため、保持する拘束のどの組み合わせについても、全ての拘束を満たす解が見つカゝるまで線形独立の小行列を計算する。あるいは、元の行列 Aが悪状態又は不十分なランクであるために完全な解を求めることが出来ない場合もあり、そのときは最も拘束を満足する解が見つ力るまで全ての小行列について計算を行う。

[0094] このように本発明では、外力及び加速度の関係が一次従属の関係を含む場合であつても、一次従属の関係を排除して線形独立の関係のみを残すことができ、物体に印加されている外力及び加速度を解析的に速やかに求めることができる。

Claims

請求の範囲

[1] 表面上の複数の点に外力が印加された場合の物体の運動を解析する運動解析方法において、

各点に印加された外力と、外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式ィ匕し、行列分解を用いて前記係数行列から非特異小行列を抽出し、抽出した非特異小行列を用いて前記連立一次方程式の解を算出し、前記連立一次方程式の途中結果及び最終結果をメモリに記録することを特徴とする運動解析方法。

[2] 前記行列分解は、 QR分解であることを特徴とする請求項 1に記載の運動解析方法

[3] QR分解を用いて前記係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて 1又は複数回の列交換及びギブンス回転を行うことを特徴とする請求項 1又は請求項 2に記載の運動解析方法。

[4] 前記上三角行列から零である対角要素を探索し、零である対角要素が探索された場合、該対角要素の右隣から順に非零である要素を探索し、非零である要素が探索された場合、該要素を含む列ベクトルと前記対角要素を含む列ベクトルとの間の列交換を行うことを特徴とする請求項 3に記載の運動解析方法。

[5] 前記物体は、床面を歩行する歩行ロボットであり、前記外力は前記床面から受ける抗カであることを特徴とする請求項 1乃至請求項 4に記載の運動解析方法。

[6] 表面上の複数の点に外力が印加された場合の物体の運動を解析する運動解析装 ¾【こ; i l /、て、

各点に印加された外力と外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式ィヒする手段と、行列分解を用いて前記係数行列から非特異小行列を抽出する手段と、抽出した非特異小行列を用いて前記連立一次方程式の解を算出する手段と、前記連立一次方程式の途中結果及び最終結果をメモリに記録する手段とを備えることを特徴とする運動解析装置。

[7] コンピュータに、表面上の複数の点に外力が印加された場合の物体の運動を解析させるコンピュータプログラムにおいて、コンピュータに、各点に印加された外力と外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式ィ匕させるステップと、コンピュータに、行列分解を用いて前記係数行列から非特異小行列を抽出させるステツプと、コンピュータに、抽出させた非特異小行列を用いて前記連立方程式の解を算出させるステップとを有することを特徴とするコンピュータプログラム。

コンピュータに、表面上の複数の点に外力が印加された場合の物体の運動を解析させるコンピュータプログラムが記録されているコンピュータでの読み取りが可能な記録媒体において、

コンピュータに、各点に印加された外力と外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式ィ匕させるステップと、コンピュータに、行列分解を用いて前記係数行列から非特異小行列を抽出させるステツプと、コンピュータに、抽出させた非特異小行列を用いて前記連立方程式の解を算出させるステップとを有するコンピュータプログラムが記録されていることを特徴とするコンピュータでの読み取りが可能な記録媒体。