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

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

Info

Publication number
JPWO2007094451A1
JPWO2007094451A1 JP2008500562A JP2008500562A JPWO2007094451A1 JP WO2007094451 A1 JPWO2007094451 A1 JP WO2007094451A1 JP 2008500562 A JP2008500562 A JP 2008500562A JP 2008500562 A JP2008500562 A JP 2008500562A JP WO2007094451 A1 JPWO2007094451 A1 JP WO2007094451A1
Authority
JP
Japan
Prior art keywords
matrix
external force
decomposition
column
motion analysis
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.)
Pending
Application number
JP2008500562A
Other languages
English (en)
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.)
Japan Science and Technology Agency
National Institute of Japan Science and Technology Agency
Original Assignee
Japan Science and Technology Agency
National Institute of Japan Science and Technology Agency
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 Japan Science and Technology Agency, National Institute of Japan Science and Technology Agency filed Critical Japan Science and Technology Agency
Publication of JPWO2007094451A1 publication Critical patent/JPWO2007094451A1/ja
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/11Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/12Simultaneous equations, e.g. systems of linear equations

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Operations Research (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • Complex Calculations (AREA)

Abstract

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

Description

本発明は、行列方程式を用いて運動解析を行う運動解析方法、運動解析装置、コンピュータプログラム、及び記録媒体に関する。
ヒューマノイドロボットの挙動を解析する動的シュミレーションでは、外部から印加される外力(接触力および衝突)を注意深く取り扱う必要がある。近年、正確かつ有効な解を見出すために精力的な研究がなされており、ペナルティ手法、撃力ベースの手法、解析的手法などの手法により運動の解析が行われている。
ペナルティ手法は、拘束を違反した度合いに応じて力を加えることで、シミュレーションを進めるうちに拘束違反を解消する手法である(例えば、非特許文献1及び2参照)。接触力を正確に求めたい場合などに用いられてきた。拘束力は、侵入量(すなわち、拘束違反の量)と速度とから直接的に求まるが、拘束が解消するまでシミュレーションを繰り返す必要があるという問題点を有している。
撃力ベースの手法では、瞬間的な撃力が物体の速度を変化させるとして衝突および接触を解析する手法である(例えば、非特許文献3参照)。この手法では、物体同士の相互作用やエネルギーの伝搬に着目しており、シュミレーションの正確さよりも、むしろ物体の挙動を視覚化することに主眼を置いている。
一方、解析的手法は、幾つかの制約があるものの現実のロボットの挙動を正確にシミュレートするのに最も適している。解析的手法では、線形拘束条件を取り扱い、拘束条件の式と運動方程式とを連立させて解くことで、拘束条件を満たすために必要な力を計算する(例えば、特許文献4及び5参照)。このように解析的手法では、拘束条件を満たす力を解くため、シミュレーションの時間をある程度広くとっても正確なシミュレーションを行うことができる。
ウィルヘルム(J. Wilhelms)、他2名著、「動的アニメーション:相互作用と制御(Dynamic animation: interaction and control)」、ビジュアルコンピュータ(Visual Computer)、1988年、第4巻、P283−295 パンネ(M. van de Panne)、他1名著、「センサ−アクチュエータ ネットワーク(Sensor-actuator network)」、プロシーディングス エーシーエム シーグラフ(Proc. ACM SIGGRAPH)、1993年、第27(2)巻、P335−343 シュミドル(H. Schmidl)、他1名著、「剛体シュミレーションのための瞬間的な接触(A fast impulsive contact suite for rigid body simulation)」、視覚化とコンピュータグラフィックス(Visualization and Computer Graphics)、2004年、第10巻、P189−197 フェザーストーン(R. Featherstone)著、「ロボット力学のアルゴリズム(Robot Dynamics Algorithms)」、クルワーアカデミック(Kluwer Academic)、1987年 バラフ(D. Baraff)著、「非侵入剛体に対する接触力の計算(Fast contact force computation for nonpenetratingrigid bodies)」、シーグラフ 94 コンファレンス(SIGGRAPH 94 Conference)、1994年、P23−34
しかしながら、前述した解析的手法では、全ての接触点について剛体の運動を表す行列方程式を解く必要があるため、接触点が多い場合などに多くの計算時間を必要とすることが問題となっていた。
また、剛体上の複数の点に作用する接触力と加速度との関係が一次独立となり、係数行列が正則でありさえすれば、LU分解法、QR分解法等の行列分解を用いることにより解を求めることが可能であるが、実際には、多くの接触状態において一次従属の拘束を含むことが知られており、係数行列について独立性を仮定することはできない。この場合、ある拘束条件を無視し、他のすべての拘束条件を満たす解を求める手法がよく採用されるが、接触力自体を求めるのと同様に複雑な計算を要するという問題点を有している。
本発明は斯かる事情に鑑みてなされたものであり、接触点が多く、かつ接触力と加速度との関係が一次従属となる場合であっても、各接触点に印加される外力及び加速度を速やかに求めることができる運動解析方法、運動解析装置、コンピュータプログラム、及び記録媒体を提供することを目的とする。
第1発明に係る運動解析方法は、表面上の複数の点に外力が印加された場合の物体の運動を解析する運動解析方法において、各点に印加された外力と、外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて前記係数行列から非特異小行列を抽出し、抽出した非特異小行列を用いて前記連立一次方程式の解を算出し、前記連立方程式の途中結果及び最終結果をメモリに記録することを特徴とする。
本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしているため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。
第2発明に係る運動解析方法は、前記行列分解は、QR分解であることを特徴とする。
本発明にあっては、行列分解としてQR分解を用いるため、既知の演算によって非特異小行列が抽出される。
第3発明に係る運動解析方法は、QR分解を用いて前記係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて1又は複数回の列交換及びギブンス回転を行うことを特徴とする。
本発明にあっては、QR分解により係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて列交換及びギブンス回転を行うようにしているため、計算コストが低減される。
第4発明に係る運動解析方法は、前記上三角行列から零である対角要素を探索し、零である対角要素が探索された場合、該対角要素の右隣から順に非零である要素を探索し、非零である要素が探索された場合、該要素を含む列ベクトルと前記対角要素を含む列ベクトルとの間の列交換を行うことを特徴とする。
本発明にあっては、上三角行列から零である対角要素と、その対角要素の右側の列ベクトルに含まれる非零の要素とを探索する。そして、両要素が探索された場合、探索された対角要素及び非零である要素を含む2つの列ベクトルを交換するようにしているため、非特異小行列が速やかに抽出される。
第5発明に係る運動解析方法は、前記物体は、床面を歩行する歩行ロボットであり、前記外力は前記床面から受ける抗力であることを特徴とする。
本発明にあっては、解析対象として歩行ロボットを採用し、歩行ロボットが床面から受ける抗力を接触力問題として解くことができる。
第6発明に係る運動解析装置は、表面上の複数の点に外力が印加された場合の物体の運動を解析する運動解析装置において、各点に印加された外力と外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化する手段と、行列分解を用いて前記係数行列から非特異小行列を抽出する手段と、抽出した非特異小行列を用いて前記連立一次方程式の解を算出する手段と、前記連立一次方程式の途中結果及び最終結果をメモリに記録する手段とを備えることを特徴とする。
本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしているため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。
第7発明に係るコンピュータプログラムは、コンピュータに、表面上の複数の点に外力が印加された場合の物体の運動を解析させるコンピュータプログラムにおいて、コンピュータに、各点に印加された外力と、外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化させるステップと、コンピュータに、行列分解を用いて前記係数行列から非特異小行列を抽出させるステップと、コンピュータに、抽出させた非特異小行列を用いて前記連立方程式の解を算出させるステップとを有することを特徴とする。
本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしているため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。
第8発明に係る記録媒体は、コンピュータに、表面上の複数の点に外力が印加された場合の物体の運動を解析させるコンピュータプログラムが記録されているコンピュータでの読み取りが可能な記録媒体において、コンピュータに、各点に印加された外力と、外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化させるステップと、コンピュータに、行列分解を用いて前記係数行列から非特異小行列を抽出させるステップと、コンピュータに、抽出させた非特異小行列を用いて前記連立方程式の解を算出させるステップとを有するコンピュータプログラムが記録されていることを特徴とする。
本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしているため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。
第1発明及び第6発明による場合は、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしている。したがって、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除して線形独立の関係のみを残すことができ、物体に印加されている外力及び加速度を解析的に速やかに求めることができる。
第2発明による場合は、行列分解として特定のQR分解を用いるため、既知の演算によって非特異小行列を抽出することができる。
第3発明による場合は、QR分解により係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて列交換及びギブンス回転を行うようにしているため、計算コストを低減することができる。
第4発明による場合は、上三角行列から零である対角要素と、その対角要素の右側の列ベクトルに含まれる非零の要素とを探索する。そして、両要素が探索された場合、探索された対角要素及び非零である要素を含む2つの列ベクトルを交換するようにしているため、非特異小行列を速やかに抽出することができる。
第5発明による場合は、解析対象として歩行ロボットを採用し、歩行ロボットが床面から受ける抗力を接触力問題として解くことができる。
第7発明及び第8発明による場合は、上述した発明に記載の運動解析方法をコンピュータプログラム及び記録媒体に記録されたコンピュータプログラムにより実現できるため、適宜のコンピュータを用いて物体に印加されている外力及び加速度を解析的に速やかに求めることができる。
本実施の形態に係る運動解析装置の内部構成を示すブロック図である。 2つの手法の計算コストの比較結果を示すグラフである。 列交換プロセスを用いた場合の計算コストを示すグラフである。 消去及び挿入プロセスを用いた場合の計算コストを示すグラフである。 非零要素の探索手法を模式的に説明する説明図である。 本発明に係るコンピュータプログラムにより実行される処理の手順を説明するフローチャートである。
符号の説明
100 情報処理装置
101 CPU
102 ROM
103 RAM
104 記憶装置
105 入出力IF
106 キーボード
107 モニタ
108 補助記憶装置
109 バス
110 記録媒体
図1は本実施の形態に係る運動解析装置の内部構成を示すブロック図である。本実施の形態に係る運動解析装置は、パーソナルコンピュータ、ワークステーション等の情報処理装置100により実現される。情報処理装置100は、演算手段としてのCPU101を備えており、このCPU101には、ROM102、RAM103、記憶装置104、入出力IF105、補助記憶装置108等のハードウェアがバス109を介して接続されている。
ROM102には、バス109に接続された各種ハードウェアの動作を制御するための制御プログラムが格納されている。CPU101は、この制御プログラムをRAM103上にロードして実行することにより、ハードウェア全体の動作を制御する。また、RAM103には、運動解析の途中結果及び最終結果が適宜書き出される。
記憶装置104は、ハードディスクドライブを備えており、本実施の形態で説明する運動解析方法を具体的に実現するためのコンピュータプログラム、及びこのコンピュータプログラムを実行する際に必要となるデータが記憶される。
入出力IF105には、入力デバイスであるキーボード106、及び出力デバイスであるモニタ107が接続される。情報処理装置100は、キーボード106を通じて運動解析の開始指示、運動解析に必要なパラメータの入力を受付ける。また、情報処理装置100は、キーボード106を通じて入力されたパラメータ、前述したコンピュータプログラムの演算結果等をモニタ107上に表示する。
補助記憶装置108は、本発明に係るコンピュータプログラムが記録されたFD、CD−ROM、DVD−ROM等の記録媒体110からコンピュータプログラムを読み取るためのFDドライブ、CD−ROMドライブ、DVD−ROMドライブ等を備えている。補助記憶装置108により読み取られたコンピュータプログラムは記憶装置104に記憶される。CPU101は、必要に応じて記憶装置104から前述のコンピュータプログラムをRAM103上にロードして実行することにより、情報処理装置100を本実施の形態に係る運動解析装置として動作させる。
情報処理装置100による運動解析の対象として歩行ロボットを採用し、歩行ロボットが床面を歩行する際の運動を解析する手法を以下で説明する。歩行ロボットは、床面からの抗力を受けるため、この抗力(接触力)と加速度との関係を運動方程式により表す。典型的な接触力の問題は、次の式のように表すことができる。
Figure 2007094451
ここで、fは物体表面の1又は複数の接触点における接触力を表すベクトルであり、aは各接触点での法線方向及び接戦方向の加速度を表すベクトルである。Aは接触力と加速度とを関連付ける行列であり、物体の質量に関係する量を表す。また、a0 は物体に接触力が印加されていない場合の前記接触点での加速度を表すベクトルである。
行列Aおよびベクトルa0 の各要素については、例えば、SDFASTのような動力学のソフトウェアを用いることにより、解析的に又は実験的に導出することができる。すなわち、行列Aおよびベクトルa0 の各要素の値は既知の量と考えることができるため、数式1は、行列Aを係数行列とした加速度aおよび接触力fに関する連立一次方程式と見ることができ、この連立一次方程式を解く問題に帰着することができる。
係数行列を用いて連立一次方程式を解く方法としては、ガウスの消去法、LU分解法、QR分解法等が知られている。このうち、QR分解に基づく方法は、係数行列が正則でありさえすれば解を求めることが可能である。しかし、行列が特異に近い場合には、解の精度を高めようとすると軸列の選択が必要なことが知られており、軸選択を行わない並列解法には適用範囲に限界がある。
実際、数式1に示した運動方程式では、多くの接触状態において一次従属の拘束を含むことが知られており、各加速度について独立性を仮定することができない。この場合、ある拘束条件を無視し、全ての拘束条件に適合する解を求める手法がよく採用されているが、接触力自体を求めるのと同様に複雑な計算を要する。また、摩擦力に起因する加速度によっては、任意の加速度の任意の次元において一次従属となり得る。
そこで、本実施の形態では、行列Aから非特異小行列を構築する手法(以下、干渉性利用法という)について説明する。この干渉性利用法では、行列Aについて消去しないサブセットの任意のインデックスを受付ける。このアルゴリズムでは、ある加速度間の拘束を無視することによって一次独立の接触問題として解を求め、解を求めた後に全ての加速度間の拘束をチェックする。不要な拘束を全て不要セットに置き換えた場合、これらの不要セットのインデックスの1つは保持セットに含まれ、行列Aの異なる一次独立の小行列を算出することができる。どの拘束が必要であるか、どの拘束が他の拘束に影響を及ぼすかということについては一般的には分からないため、保持する拘束のどの組み合わせについても、全ての拘束を満たす解が見つかるまで線形独立の小行列を計算する。あるいは、元の行列Aが悪状態又は不十分なランクであるために完全な解を求めることが出来ない場合もあり、そのときは最も拘束を満足する解が見つかるまで全ての小行列について計算を行う。
特異小行列の抽出プロセスにはQR分解の理解が不可欠であるため、まず、列交換を保持するプロセスについて説明する。QR分解では、m行n列の任意の行列Aを、m行m列の直交行列Qとm行n列の上三角行列Rとに行列分解する。
Figure 2007094451
任意の行列AをQR分解することは可能であるが、直交行列Q及び上三角行列Rは一意ではない。QR分解は、元の行列Aが調整される場合、又は最初から計算をし直す努力を避けたい場合において特に有用である。行または列が追加された場合、QR分解はアップデートされ、行または列が削除された場合、QR分解をダウンデートされる。または、ランク1の行列の加算又は減算に続いて再計算が実行される。
QR分解のアップデート及びダウンデートは、行列の任意の列の交換に続いて分解を保持するという標準的な手法が採用される。例えば、第n列と第n+k列との間の交換を行う場合、第n列および第n+k列を消去し、その新たな2箇所に列を挿入し、それぞれ2回のアップデート及びダウンデートを行う。しかしながら、これらの演算は単一列の交換演算により統合することができ、演算資源をかなり削減することができる。
QR分解のアップデート及びダウンデートはギブンス回転を用いて実行される。ギブンス行列Gn(i,j,θ)は、n次元空間のうちの2つの部分空間(1つはi次元の部分空間、もう1つはj次元の部分空間とする)の間の回転を表している。ギブンス行列Gn(i,j,θ)は、n×nの単位行列のうち、(i,i)、(i,j)、(j,i)、(j,j)の各要素をそれぞれcosθ、−sinθ、sinθ、cosθにより置き換えた行列である。QR分解は以下のように与えられる。このとき、Q’の直交性は保持されるが、Rは上三角行列でなくなる。
Figure 2007094451
行列AのQR分解が与えられた場合、以下のようにして第k列を消去して行列分解を保持することができる。
Figure 2007094451
例えば、6行6列の行列Rから第3列の削除を考えた場合、行列Rの構造は以下のように影響を受ける。
Figure 2007094451
ここで、「×」のシンボル、及び「#1 」〜「#3 」のシンボルで表される各要素は一般に非零の値をとる(但し、偶然的に零の値をとり得る)。また、シンボルが付されていない各要素は零の値をとる。行列R(1) は現時点で上三角行列ではないが、ギブンス回転を3回適用することによって対角成分の直下の要素(すなわち、「#1 」〜「#3 」のシンボルで表した各要素)を1つずつ零にすることができ、行列R(1) の上三角性を修復することができる。ここで、ギブンス回転は、インデックスを付した順に隣り合う行同士について適用され、その結果、対角要素の下側に非零の要素が生成されないことが保証される。
行列Aに列を挿入し、QR分解をアップデートする場合も同様である。行列AのQR分解が与えられた場合、以下のようにして、行列Aに新しい列cT を追加し、行列分解を保持することができる。
Figure 2007094451
例えば、5行4列の新たに第3列を挿入することを考えた場合、行列Rの構造は以下のように影響を受ける。
Figure 2007094451
上の式においても「×」のシンボル、及び「#1 」〜「#3 」のシンボルで表される各要素は一般に非零の値をとる(但し、偶然的に零の値をとり得る)。また、「・」のシンボルは、元々は零の要素であるが、ギブンス回転を適用することによって非零の値に変化し得る要素を示している。例えば、「#1」の要素を零にするためにギブンス回転を適用した場合、「・1」の要素が非零の値になり、「#2」の要素を零にするためにギブンス回転を適用した場合、「・2」の要素が非零の値になる。したがって、この場合には行列R(1) を上三角行列に変換するために2回のギブンス回転を行う必要がある。
同様の手法は、行列Rから2つの列を入れ替えた結果として生成される行列R(1) の構造を修復する際にも適用される。2つの列の入れ替えは、列の挿入と削除とを同時的に実行することと等価であり、R(1) の構造を修復する前にアップデートと組み合わせて複数回の列の消去と挿入とを実行する場合と比べた場合、演算資源を削減する効果が高い。行列Aの第i列と第j列との入れ替えを考えた場合、以下のように行列分解を保持することができる。
Figure 2007094451
行列R(1) の構造はもはや上三角行列ではない。8行8列の上三角行列を例にとり、その第3列と第6列とを入れ替えた行列R(1) の構造を修復するプロセスについて以下に説明する。行列R(1) は以下のように表すことができる。
Figure 2007094451
この例における行列R(1) は、上三角行列の第3列と第6列とを入れ替えた行列であるから、その第3列には対角要素より下側に3つの非零の要素が存在する。これらの要素を零にするためにギブンス回転を行う。このとき、ギブンス回転の影響を受ける要素の数を最小にするために、隣り合う行同士にギブンス回転を適用する。「#1 」のシンボルで示した要素を零にするためにギブンス回転を行った結果、「・1 」のシンボルで示した要素は零から非零の値に変わる。「・2 」、「・3 」のシンボルで示した要素についても同様である。
変形した結果の行列をR(2) とした場合、行列R(2) は以下のように表すことができる。
Figure 2007094451
上の式で示したように、行列R(2) には、第4列及び第5列の対角要素の下側に「#4 」及び「#5 」のシンボルで示した非零の要素が存在する。この第2段階での変形では、第1段階の3回のギブンス回転により生成されたこれらの非零の要素を零にするために、2回のギブンス回転を行う。その結果、「#4 」及び「#5 」で示した要素は零となり、「・4 」及び「・5 」で示した要素はそのときのギブンス回転の影響を受けて値が変化する。この第2段階での変形を行った結果、前述の行列R(2) は上三角行列R(3) に変形される。対応する直交行列Q(3) は、直交行列Qに逆の回転演算を施すことにより得られる。
次に、計算コストについて説明する。
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)−1となる。
一方、列交換の手法を用いる場合、前述したように、第1段階でk回のギブンス回転を行い、第2段階でk−1回のギブンス回転を行うため、修復するために要するギブンス回転の数は2k−1でよい。すなわち、トータルで4(N−(n+k))回のギブンス回転を節約することができる。交換すべき2つの列のインデックスのうち、大きい方のインデックスにより節約できるギブンス回転の数が定まり、インデックスが元の行列の最後の列に近づくにつれて減少する。
しかしながら、ギブンス回転は非零の要素を含んだ列に対してのみ適用すればよいので、前述した計算は必ずしも真の節約量を反映していない。この点を考慮し、行列回転を用いるよりも列回転を用いて節約量を計算することができる。N行N列の行列について第n列と第n+k列との間の列交換を行う場合、以下の回数だけギブンス列回転が必要となる。
Figure 2007094451
同様に、M行N列の行列から列の消去を行う場合、以下の回数だけギブンス列回転が必要となる。但し、Mは、M∈[N,N+1]を満足する整数である。
Figure 2007094451
更に、列の挿入を行う場合には、以下の回数だけギブンス列回転が必要となる。但し、Mは、M∈[N+1,N+2]を満足する整数である。
Figure 2007094451
したがって、2つの列を消去して2つの列を挿入する場合の計算コストは、以下のように表すことができる。
Figure 2007094451
図2は2つの手法の計算コストの比較結果を示すグラフである。このグラフでは24行24列の行列について計算コストを算出しており、底面の2つの軸はそれぞれ交換対象の列のインデックスを示し、高さ方向の軸は相対コストを示している。このグラフから、消去及び挿入プロセスに比べて列交換プロセスは大部分の条件下で計算コストが良いことが分かる。相対コストが1となる条件をグラフ中の破線により示している。すなわち、交換対象の列のインデックスの何れか一方が20以下である場合、計算コストは交換プロセスを使用する方が良く、交換対象の列のインデックスが元の行列の最後の列(第24列)に近づくにつれて計算コストは消去及び挿入プロセスを使用する方が良くなる。
図3は列交換プロセスを用いた場合の計算コストを示すグラフであり、図4は消去及び挿入プロセスを用いた場合の計算コストを示すグラフである。何れのグラフについても24行24列の行列について計算コストを算出している。何れのプロセスにおいても計算コストは元の行列のサイズに依存する。行列サイズが大きくなるにつれ、列交換プロセスの計算コストは、消去及び挿入プロセスの計算コストの50%に近づくことが分かった。このファクターは、任意の2つの列の交換の可能性を等確率と仮定し、列回転を実行する回数をカウントすることにより見積もった値である。4行4列の行列、又はそれ以上のサイズを有する行列では、列交換プロセスは消去及び挿入プロセスよりも計算速度が平均的に大きく、3行3列の行列で計算コストが等しくなる。
次に、最大の非特異小行列を抽出する手法について説明する。
与えられた行列から最大限に大きな非特異小行列を抽出する手法は、前述したようなQR分解を基にしており、以下の仮定を必要とする。まず、与えられた行列が対象行列、又は行交換した対象行列であることを仮定する。また、正規直交のQR分解において行列Rが与えられていることを仮定する。更に、消去することができない行および列の組を示す保持用インデックスのセットが与えられていることを仮定する。更に、保持用インデックスのセットの軽微な変更(典型的には、新たなインデックスの追加)に続いて小行列の抽出を行うことを仮定する。更に、一連の類似の行列について一から小行列の抽出を繰り返すことを仮定する。
行列の左上側のブロックの非特異小行列と明確に区別するために元の行列の行および列を再配列することを行う。再配列するインデックスは、続くプロセスの繰り返しの間だけ内部的に保持される。
前述の仮定を用いて非特異小行列を抽出するために、対象行列Aに対するQR分解を以下のように表す。
Figure 2007094451
ここで、0又は3のインデックスを付した小行列(例えば、A0 、A3 等 )は正方行列であり、同じインデックスを付した全ての小行列は同じ次元を有する。小行列R0 の対角要素の全てを非零と仮定したとき、小行列A0 が非特異であり、かつ、小行列A0 の外部に行または列の組を追加することによって拡張した小行列A0 は特異行列となることが証明できる。
QR分解を構築することは、行列Aの各列に亘って互いに直交するベクトルのセットを見つけることと等価である。また、QR分解を区分化することは、基を構成するまで直交するベクトルのセットを縮小することと等価である。すなわち、行列Rの左上隅から完全なシリーズとなるように、行列Rの対角要素にできるだけ多くの非零の要素を配列させるプロセスが要求される。
そのため、行列Rの対角要素を左上隅から順に調べ、零の要素を探索する。本実施の形態では、微小な値を閾値として予め設定しておき、その閾値より小さい値を見つけた場合、その要素は零であると判定する。零の要素を見つけた場合、その要素のすぐ右隣から列のインデックスが増える方向に非零の要素を探索する。その行に非零の要素が見つからない場合には、行を一段だけ下げて対角要素の右側を探索し、非零の要素が見つかるまで行を一段ずつ下げながら行列Rの右下隅に至るまで探索を繰り返す。非零要素が見つからない場合、行列Rは区分化された構造を有しているため、本プロセスを終了する。なお、対角要素の下側は零であることが保証されているため、探索を省略できることはいうまでもない。図5は非零要素の探索手法を模式的に説明する説明図である。
このような要素が見つかった場合、その要素を含む列は零対角要素を含んだ列と交換され、アップデートされた行列Rの次の対角要素が零であるか否かの判定が行われる。列交換は、非零対角要素の上側及び左側には影響を及ぼさないため、行列Rの構造は保持される。小行列における零対角要素の下側及び右側に非零要素が存在するということは、その零対角要素を含んだ列と行列Rの非零対角要素を含んだ各列とが一次従属でないことを保証する。
保持用セットは、行列Rの構造を修復するのに先立ち、保持した列を非特異ブロックの最初の列と交換することにより適用することができる。例えば、ある列が保持用セットに存在するが非特異ブロックに存在しない場合、その列はそれ自身が保持用セットに含まれない非特異ブロックの最初の列と交換される。非特異ブロック中の目的の列の対角要素が零になるという結果である場合には、保持セットは一次従属の関係にある列を含み、もはや解が存在しないためプロセスを終了する。
この列交換の手続きは、保持用セットの操作は最小の列交換アップデートにより処理することができ、アルゴリズムの次の反復の間で期待された干渉性を使用することができることを確実にする。また、保持用セットが実現可能である場合には速やかに決定することができる。このパターンを用いてオリジナルのマトリックスの順序を入れ替えることにより、最初のQR分解は区分化した形態になるか、又は近づき、必要な列交換の数を減らすことができる。
このようにして求めた解は対象行列と同等に扱うことができるが、結果として導出された非正則ブロックはアルゴリズムの実行を繰り返す度にチェックする必要がある。これは、候補となるブロック行列を完全にQR分解することを要求する。Q行列がアルゴリズムを通して保持される場合、分解を再計算する必要はなくなるが、候補となる行列のオリジナルの行列に対する相対サイズに依存することとなり、Q行列の間接的な保持及びQR分解を一から再計算することを避けることが十分可能になる。
次に、本発明に係るコンピュータプログラムにより実行される処理の手順を説明する。図6は本発明に係るコンピュータプログラムにより実行される処理の手順を説明するフローチャートである。外部から入力されたパラメータに基づいて物体(例えば、歩行ロボット)の運動解析を行う場合、物体表面上の各点に印加された外力と加速度との関係を定式化する(ステップS1)。物体表面上の各点に印加された外力と加速度との関係は、前述の数式1に示したように連立一次方程式を用いることにより、物体表面の1又は複数の接触点における接触力を表すベクトルf、各接触点での法線方向及び接線方向の加速度を表すベクトルa、接触力と加速度とを関連付ける行列Aにより記述することができる。
次いで、行列分解(特に、QR分解)を用いることにより行列Aから非特異小行列を抽出し(ステップS2)、抽出した非特異小行列を用いて定式化した連立一次方程式の解を算出する(ステップS3)。このとき、ある加速度間の拘束を無視することによって一次独立の接触問題として解を求め、解を求めた後に全ての加速度間の拘束をチェックする。不要な拘束を全て不要セットに置き換えた場合、これらの不要セットのインデックスの1つは保持セットに含まれ、行列Aの異なる一次独立の小行列を算出することができる。どの拘束が必要であるか、どの拘束が他の拘束に影響を及ぼすかということについては一般的には分からないため、保持する拘束のどの組み合わせについても、全ての拘束を満たす解が見つかるまで線形独立の小行列を計算する。あるいは、元の行列Aが悪状態又は不十分なランクであるために完全な解を求めることが出来ない場合もあり、そのときは最も拘束を満足する解が見つかるまで全ての小行列について計算を行う。
このように本発明では、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除して線形独立の関係のみを残すことができ、物体に印加されている外力及び加速度を解析的に速やかに求めることができる。
【0003】
グラム、及び記録媒体を提供することを目的とする。
課題を解決するための手段
[0009]
第1発明に係る運動解析方法は、表面上の複数の点に外力が印加された場合の物体の運動を解析する運動解析方法において、各点に印加された外力と、外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化し、QR分解を用いて前記係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて1又は複数回の列交換及びギブンス回転を行うことにより、前記係数行列から非特異小行列を抽出し、抽出した非特異小行列を用いて前記連立一次方程式の解を算出し、前記連立一次方程式の途中結果及び最終結果をメモリに記録することを特徴とする。
[0010]
本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしているため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。また、QR分解により係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて列交換及びギブンス回転を行うようにしているため、計算コストが低減される。
[0011]
第2発明に係る運動解析方法は、表面上の複数の点に外力が印加された場合の物体の運動を解析する運動解析方法において、各点に印加された外力と、外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化し、QR分解を用いて前記係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて列交換を行う処理と列交換によって対角要素より下側に非零の要素が存在することとなった場合、前記非零の要素を含んだ行ベクトル及び該行ベクトルに隣接する行ベクトルについてギブンス回転を行う処理とを繰り返すことにより、前記係数行列から非特異小行列を抽出し、抽出した非特異小行列を用いて前記連立一次方程式の解を算出し、前記連立一次方程式の途中結果及び最終結果をメモリに記録することを特徴とする。
[0012]
本発明にあっては、QR分解を用いて係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて列交換を行う処理と列交換によって対角要素より下側に非零の要素が存在することとなった場合、非零の要素を含んだ行ベクトル及びこの行ベクトルに隣接する行ベクトルについてギブンス回転を行う処理とを繰り返すようにしているため、計算コストが低減される。
[0013]
[0014]
[0015]
第4発明に係る運動解析方法は、前記上三角行列から零である対角要素を探索し、零である対角要素が探索された場合、該対角要素の右隣から順に非零である要素を探索し、非零である要素が探索された場合、該要素を含む列ベクトルと前記対角要素を含む列ベクトルとの間の列交換を行うことを特徴とする。
[0016]
本発明にあっては、上三角行列から零である対角要素と、その対角要素の右側の列ベクトルに含まれる非零の要素とを探索する。そして、両要素が探索された場合、
【0004】
探索された対角要素及び非零である要素を含む2つの列ベクトルを交換するようにしているため、非特異小行列が速やかに抽出される。
[0017]
第5発明に係る運動解析方法は、前記物体は、床面を歩行する歩行ロボットであり、前記外力は前記床面から受ける抗力であることを特徴とする。
[0018]
本発明にあっては、解析対象として歩行ロボットを採用し、歩行ロボットが床面から受ける抗力を接触力問題として解くことができる。
[0019]
第6発明に係る運動解析装置は、上述した発明の何れか1つに記載の運動解析方法を用いて物体の運動を解析することを特徴とする。
[0020]
本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしているため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。また、QR分解により係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて列交換及びギブンス回転を行うようにしているため、計算コストが低減される。
[0021]
第7発明に係るコンピュータプログラムは、上述した発明の何れか1つに記載の運動解析方法を用いて物体の運動を解析させるステップを有することを特徴とする。
[0022]
本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしているため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。また、QR分解により係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて列交換及びギブンス回転を行うようにしているため、計算コストが低減される。
【0005】
[0023]
第8発明に係る記録媒体は、上述した発明の何れか1つに記載の運動解析方法を用いて物体の運動を解析させるステップを有するコンピュータプログラムが記録されていることを特徴とする。
[0024]
本発明にあっては、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしているため、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除することにより線形独立の関係のみが残る。また、QR分解により係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて列交換及びギブンス回転を行うようにしているため、計算コストが低減される。
発明の効果
[0025]
第1発明及び第6発明による場合は、外力と加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて係数行列から非特異小行列を抽出するようにしている。したがって、外力及び加速度の関係が一次従属の関係を含む場合であっても、一次従属の関係を排除して線形独立の関係のみを残すことができ、物体に印加されている外力及び加速度を解析的に速やかに求めることができる。
[0026]
また、本発明による場合は、行列分解として特定のQR分解を用いるため、既知の演算によって非特異小行列を抽出することができる。
[0027]
第2発明による場合は、QR分解を用いて係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて列交換を行う処理と列交換によって対角要素より下側に非零の要素が存在することとなった場合、非零の要素を含んだ行ベクトル及びこの行ベクトルに隣接する行ベクトルについてギブンス回転を行う処理とを繰り返すようにしているため、計算コストを低減することができる。
[0028]
第4発明による場合は、上三角行列から零である対角要素と、その対角要素の右側の列ベクトルに含まれる非零の要素とを探索する。そして、両要素が探索された場合、探索された対角要素及び非零である要素を含む2つの列ベクトルを交換するようにしているため、非特異小行列を速やかに抽出することができる。

Claims (8)

  1. 表面上の複数の点に外力が印加された場合の物体の運動を解析する運動解析方法において、
    各点に印加された外力と、外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化し、行列分解を用いて前記係数行列から非特異小行列を抽出し、抽出した非特異小行列を用いて前記連立一次方程式の解を算出し、前記連立一次方程式の途中結果及び最終結果をメモリに記録することを特徴とする運動解析方法。
  2. 前記行列分解は、QR分解であることを特徴とする請求項1に記載の運動解析方法。
  3. QR分解を用いて前記係数行列を直交行列と上三角行列との積に分解し、得られた上三角行列を構成する列ベクトルについて1又は複数回の列交換及びギブンス回転を行うことを特徴とする請求項1又は請求項2に記載の運動解析方法。
  4. 前記上三角行列から零である対角要素を探索し、零である対角要素が探索された場合、該対角要素の右隣から順に非零である要素を探索し、非零である要素が探索された場合、該要素を含む列ベクトルと前記対角要素を含む列ベクトルとの間の列交換を行うことを特徴とする請求項3に記載の運動解析方法。
  5. 前記物体は、床面を歩行する歩行ロボットであり、前記外力は前記床面から受ける抗力であることを特徴とする請求項1乃至請求項4に記載の運動解析方法。
  6. 表面上の複数の点に外力が印加された場合の物体の運動を解析する運動解析装置において、
    各点に印加された外力と外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化する手段と、行列分解を用いて前記係数行列から非特異小行列を抽出する手段と、抽出した非特異小行列を用いて前記連立一次方程式の解を算出する手段と、前記連立一次方程式の途中結果及び最終結果をメモリに記録する手段とを備えることを特徴とする運動解析装置。
  7. コンピュータに、表面上の複数の点に外力が印加された場合の物体の運動を解析させるコンピュータプログラムにおいて、
    コンピュータに、各点に印加された外力と外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化させるステップと、コンピュータに、行列分解を用いて前記係数行列から非特異小行列を抽出させるステップと、コンピュータに、抽出させた非特異小行列を用いて前記連立方程式の解を算出させるステップとを有することを特徴とするコンピュータプログラム。
  8. コンピュータに、表面上の複数の点に外力が印加された場合の物体の運動を解析させるコンピュータプログラムが記録されているコンピュータでの読み取りが可能な記録媒体において、
    コンピュータに、各点に印加された外力と外力が印加された場合の各点での加速度との関係を係数行列を用いた連立一次方程式により定式化させるステップと、コンピュータに、行列分解を用いて前記係数行列から非特異小行列を抽出させるステップと、コンピュータに、抽出させた非特異小行列を用いて前記連立方程式の解を算出させるステップとを有するコンピュータプログラムが記録されていることを特徴とするコンピュータでの読み取りが可能な記録媒体。
JP2008500562A 2006-02-16 2007-02-16 運動解析方法、運動解析装置、コンピュータプログラム、及び記録媒体 Pending JPWO2007094451A1 (ja)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
JP2006039834 2006-02-16
JP2006039834 2006-02-16
PCT/JP2007/052845 WO2007094451A1 (ja) 2006-02-16 2007-02-16 運動解析方法、運動解析装置、コンピュータプログラム、及び記録媒体

Publications (1)

Publication Number Publication Date
JPWO2007094451A1 true JPWO2007094451A1 (ja) 2009-07-09

Family

ID=38371627

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2008500562A Pending JPWO2007094451A1 (ja) 2006-02-16 2007-02-16 運動解析方法、運動解析装置、コンピュータプログラム、及び記録媒体

Country Status (2)

Country Link
JP (1) JPWO2007094451A1 (ja)
WO (1) WO2007094451A1 (ja)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP6421683B2 (ja) * 2015-04-14 2018-11-14 トヨタ自動車株式会社 最適制御装置、最適制御方法及び最適制御プログラム
CN107186753B (zh) * 2017-05-17 2020-10-09 上海电器科学研究所(集团)有限公司 工业机器人性能测试的工作空间确定方法
CN107844460B (zh) * 2017-07-24 2020-12-25 哈尔滨工程大学 一种基于p-maxq的多水下机器人的围捕方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6993193B2 (en) * 2002-03-26 2006-01-31 Agilent Technologies, Inc. Method and system of object classification employing dimension reduction

Also Published As

Publication number Publication date
WO2007094451A1 (ja) 2007-08-23

Similar Documents

Publication Publication Date Title
Palar et al. On efficient global optimization via universal Kriging surrogate models
JP5250076B2 (ja) 構造予測モデル学習装置、方法、プログラム、及び記録媒体
Yuan et al. A novel adaptive importance sampling algorithm based on Markov chain and low-discrepancy sequence
JP5282493B2 (ja) 最適解関係表示装置、方法、及びプログラム
US11886783B2 (en) Simulation system for semiconductor process and simulation method thereof
US20190034802A1 (en) Dimensionality reduction in Bayesian Optimization using Stacked Autoencoders
JP2007329415A (ja) データ処理方法、データ処理プログラム、該プログラムを記録した記録媒体およびデータ処理装置
JP5402351B2 (ja) 多目的最適化設計支援装置、方法、及びプログラム
JP6381962B2 (ja) シミュレーションシステム及び方法と該システムを含むコンピュータシステム
JP7283065B2 (ja) 推定装置、最適化装置、推定方法、最適化方法、及びプログラム
Nitzler et al. A generalized probabilistic learning approach for multi-fidelity uncertainty quantification in complex physical simulations
JPWO2007094451A1 (ja) 運動解析方法、運動解析装置、コンピュータプログラム、及び記録媒体
CN113112414A (zh) 噪声估计方法、噪声估计程序以及噪声估计设备
Van Steenkiste et al. Sensitivity analysis of expensive black-box systems using metamodeling
JP2008015841A (ja) 回路解析方法、及び回路解析プログラム、回路シミュレーション装置
JP5888782B2 (ja) 連立一次方程式の演算処理システム
JP5813498B2 (ja) モデル学習装置、関連情報抽出装置、関連情報予測装置、モデル学習方法、関連情報抽出方法、関連情報予測方法およびプログラム
JP7452648B2 (ja) 学習方法、学習装置及びプログラム
Tsilifis et al. Dimensionality Reduction for Multi-Fidelity Gaussian Processes using Bayesian Adaptation
Ayubian et al. GPU-based monte-carlo simulation for a sea ice load application
JP2006119839A (ja) 線形連立一次方程式の数値解析方法を用いた理工学シミュレータ
JP4061180B2 (ja) 線形連立1次方程式の数値解析方法を用いた理工学シミュレータ
JP7420148B2 (ja) 学習装置、学習方法及びプログラム
US20220414184A1 (en) Data processing apparatus and data processing method
CN113609720B (zh) 一种有限元分析的主从自由度处理方法、设备及存储介质

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20090422

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20110705

A02 Decision of refusal

Free format text: JAPANESE INTERMEDIATE CODE: A02

Effective date: 20111108