JP2007213172A - モンテカルロ評価における危険不確定パラメータの検出方法 - Google Patents
モンテカルロ評価における危険不確定パラメータの検出方法 Download PDFInfo
- Publication number
- JP2007213172A JP2007213172A JP2006030256A JP2006030256A JP2007213172A JP 2007213172 A JP2007213172 A JP 2007213172A JP 2006030256 A JP2006030256 A JP 2006030256A JP 2006030256 A JP2006030256 A JP 2006030256A JP 2007213172 A JP2007213172 A JP 2007213172A
- Authority
- JP
- Japan
- Prior art keywords
- uncertain
- uncertain parameter
- monte carlo
- risk
- distribution
- 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
Links
Images
Landscapes
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
【課題】再度のモンテカルロ・シミュレーションを実施することなく、評価結果に重大な影響を及ぼす危険不確定パラメータを検出する方法を提供する。
【解決手段】失敗例に対応するn個の不確定パラメータ・ベクトルΕF k(k=1,…,n)を抽出し、第j不確定パラメータから構成される第j部分集合についての平均ε'F(j)および分散s2(j)を確率変数とする新たな確率分布を生成し、その平均および分散を第j母集合の平均μ(j)および分散σ2(j)に基づいて標準化し、それぞれz(j)、χ2(j)とし、z(j)分布における裾の累積面積Pz(j)とχ2(j)分布における裾の累積面積Pχ2(j)をそれぞれ算出し、これらの代表値であるP(j)=min{Pz(j),Pχ2(j)}を、設定した危険率αで検定を行う。
【選択図】図1
【解決手段】失敗例に対応するn個の不確定パラメータ・ベクトルΕF k(k=1,…,n)を抽出し、第j不確定パラメータから構成される第j部分集合についての平均ε'F(j)および分散s2(j)を確率変数とする新たな確率分布を生成し、その平均および分散を第j母集合の平均μ(j)および分散σ2(j)に基づいて標準化し、それぞれz(j)、χ2(j)とし、z(j)分布における裾の累積面積Pz(j)とχ2(j)分布における裾の累積面積Pχ2(j)をそれぞれ算出し、これらの代表値であるP(j)=min{Pz(j),Pχ2(j)}を、設定した危険率αで検定を行う。
【選択図】図1
Description
本発明は、モンテカルロ評価における危険不確定パラメータの検出方法、特に、再度のモンテカルロ評価を実施することなくシステムの性能に重大な影響を及ぼす危険不確定パラメータを短時間で検出しモンテカルロ評価における計算負荷を大幅に軽減することが可能なモンテカルロ評価における危険不確定パラメータの検出方法に関するものである。
モンテカルロ・シミュレーションは、実環境での試験が容易ではなく、かつ、環境条件・システムの特性・初期条件等様々な不確定パラメータの存在下で動作することが要求されるシステムを事前に評価するための、乱数を利用した数値シミュレーション技術である(例えば、図14を参照。)。近年の計算機能力の大幅な向上と相俟って、モンテカルロ・シミュレーションは、航空宇宙機、船舶、各種プラント等の様々な分野におけるシステムの評価に適している。
ところで、システムの評価結果が望ましくない場合には、その原因を特定して改善処置を行う必要がある。しかし、モンテカルロ・シミュレーションによる評価では、多数の不確定パラメータを同時に加えるため、評価結果に主たる影響を及ぼす危険不確定パラメータを検出することは容易ではない。特に、2つ以上のパラメータが複合的に作用してシステムに悪影響を及ぼしている場合は、更に困難となる。
下記の非特許文献1および2において提案される方法は、上記技術的課題を解決し、不確定パラメータが複合的に作用している場合であっても、問題となる危険不確定パラメータを見つけ出すことができる。しかし、これらの方法では、問題となる危険不確定パラメータを見つけ出すために、不確定パラメータ・ベクトルを加工した後に再度モンテカルロ・シミュレーション評価を実施する必要がある。ところが、モンテカルロ・シミュレーションは計算負荷が大きい技法であり、評価対象となるシステムによっては数日から数週間の計算時間を要し、時間的コストが大きい。
ところで、システムの評価結果が望ましくない場合には、その原因を特定して改善処置を行う必要がある。しかし、モンテカルロ・シミュレーションによる評価では、多数の不確定パラメータを同時に加えるため、評価結果に主たる影響を及ぼす危険不確定パラメータを検出することは容易ではない。特に、2つ以上のパラメータが複合的に作用してシステムに悪影響を及ぼしている場合は、更に困難となる。
下記の非特許文献1および2において提案される方法は、上記技術的課題を解決し、不確定パラメータが複合的に作用している場合であっても、問題となる危険不確定パラメータを見つけ出すことができる。しかし、これらの方法では、問題となる危険不確定パラメータを見つけ出すために、不確定パラメータ・ベクトルを加工した後に再度モンテカルロ・シミュレーション評価を実施する必要がある。ところが、モンテカルロ・シミュレーションは計算負荷が大きい技法であり、評価対象となるシステムによっては数日から数週間の計算時間を要し、時間的コストが大きい。
Motoda,T.and Miyazawa Y.,"Identification of Influential Uncertainties in Monte Carlo Analysis," Journal of Spacecraft and Rockets,V0L39,No.4,2002,pp.615-623.
Motoda,T.,"Simplified Approach to Identifying Influential Uncertainties in Monte Carlo Analysis,"Journal of Spacecraft and Rockets,V0L39,No.4,2004,pp.1071-1075.
上述した通り、上記従来のモンテカルロ・シミュレーションにおける危険不確定パラメータを検出する方法は、不確定パラメータが複合的に作用している場合でも、問題となる不確定パラメータを見つけ出すことができる。
しかしながら、これらの方法では、目的とする不確定パラメータを見つけ出すために、不確定パラメータに対し適当な処置を施した後に再度モンテカルロ・シミュレーションを実施する必要があるため、時間的コストが大きく、コスト及び作業効率の面からも改善が望まれている。
そこで、本発明は、上記実情に鑑み創案されたものであって、再度のモンテカルロ評価を実施することなくシステムの性能に重大な影響を及ぼす危険不確定パラメータを短時間で検出しモンテカルロ評価における計算負荷を大幅に軽減することが可能なモンテカルロ評価における危険不確定パラメータの検出方法を提供することを目的とする。
しかしながら、これらの方法では、目的とする不確定パラメータを見つけ出すために、不確定パラメータに対し適当な処置を施した後に再度モンテカルロ・シミュレーションを実施する必要があるため、時間的コストが大きく、コスト及び作業効率の面からも改善が望まれている。
そこで、本発明は、上記実情に鑑み創案されたものであって、再度のモンテカルロ評価を実施することなくシステムの性能に重大な影響を及ぼす危険不確定パラメータを短時間で検出しモンテカルロ評価における計算負荷を大幅に軽減することが可能なモンテカルロ評価における危険不確定パラメータの検出方法を提供することを目的とする。
前記目的を達成するために、請求項1に記載のモンテカルロ評価における危険不確定パラメータの検出方法は、モンテカルロ法を用いたシミュレーションにおいてその評価結果に重大な影響を及ぼす危険不確定パラメータを検出する方法であって、乱数によって発生させるM個の不確定パラメータの内で第j不確定パラメータをε(j)(j=1,…,M)と、全シミュレーションの回数をNと、第k回目のシミュレーションにおいて発生させたM個の不確定パラメータの内で第j不確定パラメータをεk(j)(j=1,…,M)とする時、該不確定パラメータを成分に持つ不確定パラメータ・ベクトルΕk=t[εk(j)|j=1,…,M](k=1,…,N)を要素に持つ母集合[Εk|k=1,…,N]から、失敗例に対応するn個の不確定パラメータ・ベクトルΕF k(k=1,…,n)を抽出し、該抽出した各不確定パラメータ・ベクトルの第j成分εF k(j)(k=1,…,n)のみから構成される第j部分集合Ω'=[εF k(j)|k=1,…,n]の平均ε'F(j)および分散s2(j)を新たな確率変数とする確率分布を求め、前記母集合の第j不確定パラメータεk(j)(k=1,…,N)のみから構成される第j母集合Ω=[εk(j)|k=1,…,N]の確率分布に基づいて、前記第j部分集合Ω'に係る確率変数の分布に対し統計的検定を適用することにより評価結果に重大な影響を及ぼす危険不確定パラメータを求めることを特徴とする。
本願発明者が失敗ケースの不確定パラメータの特性について鋭意研究した結果、失敗ケースにおける不確定パラメータ・ベクトルの各成分から構成される部分集合の確率分布は、全ケースの不確定パラメータ・ベクトルの各成分から構成される母集合の確率分布に対して、母集合から無作為に抽出したとは考えられない程、確率分布の形態が異なることを見出した。
そこで、上記モンテカルロ評価における危険不確定パラメータの検出方法では、失敗ケースの不確定パラメータ・ベクトルの第j成分から構成される第j部分集合Ω'(j=1,…,M)を作成し、その第j部分集合Ω'についての新たな確率分布に対し、同じく第j成分から構成される第j母集合Ω(j=1,…,M)に係る確率分布に基づいて統計的検定を実施することにより、再度のモンテカルロ評価を実施することなく、システムの評価結果に重大な影響を及ぼす危険不確定パラメータを短時間で検出することが可能となる。これにより、計算時間が大幅に短縮するようになる。
本願発明者が失敗ケースの不確定パラメータの特性について鋭意研究した結果、失敗ケースにおける不確定パラメータ・ベクトルの各成分から構成される部分集合の確率分布は、全ケースの不確定パラメータ・ベクトルの各成分から構成される母集合の確率分布に対して、母集合から無作為に抽出したとは考えられない程、確率分布の形態が異なることを見出した。
そこで、上記モンテカルロ評価における危険不確定パラメータの検出方法では、失敗ケースの不確定パラメータ・ベクトルの第j成分から構成される第j部分集合Ω'(j=1,…,M)を作成し、その第j部分集合Ω'についての新たな確率分布に対し、同じく第j成分から構成される第j母集合Ω(j=1,…,M)に係る確率分布に基づいて統計的検定を実施することにより、再度のモンテカルロ評価を実施することなく、システムの評価結果に重大な影響を及ぼす危険不確定パラメータを短時間で検出することが可能となる。これにより、計算時間が大幅に短縮するようになる。
請求項2に記載のモンテカルロ評価における危険不確定パラメータの検出方法では、前記統計的検定は、前記第j部分集合Ω'の平均ε'F(j)(j=1,…,M)を各jについて求め、次いで前記第j母集合Ωの平均μ(j)(j=1,…,M)および分散σ2(j)(j=1,…,M)を用いて、
(ε'F(j)−μ(j))/(σ2(j)/n)1/2(j=1,…,M)を新たな第1確率変数z(j)(j=1,…,M)とする確率分布を生成し、該第1確率変数z(j)についての確率密度関数を用いて該第1確率変数z(j)分布の裾の第1累積面積Pz(j)(j=1,…,M)を求めることから成ることとした。
上記モンテカルロ評価における危険不確定パラメータの検出方法では、上記構成とすることにより、z≧|z(j)|になる確率Pz(j)が求まる。その確率が小さい程、現実に起こった第j部分集合Ω'の平均ε'F(j)は母集団からn個を無作為に抽出した場合には起こりにくいことを表す。つまり、確率Pz(j)の値が小さい程、それに対応する不確定パラメータε(j)は、失敗ケースに関わっている可能性が高いことになる。
(ε'F(j)−μ(j))/(σ2(j)/n)1/2(j=1,…,M)を新たな第1確率変数z(j)(j=1,…,M)とする確率分布を生成し、該第1確率変数z(j)についての確率密度関数を用いて該第1確率変数z(j)分布の裾の第1累積面積Pz(j)(j=1,…,M)を求めることから成ることとした。
上記モンテカルロ評価における危険不確定パラメータの検出方法では、上記構成とすることにより、z≧|z(j)|になる確率Pz(j)が求まる。その確率が小さい程、現実に起こった第j部分集合Ω'の平均ε'F(j)は母集団からn個を無作為に抽出した場合には起こりにくいことを表す。つまり、確率Pz(j)の値が小さい程、それに対応する不確定パラメータε(j)は、失敗ケースに関わっている可能性が高いことになる。
請求項3に記載のモンテカルロ評価における危険不確定パラメータの検出方法では、前記統計的検定は、前記第j部分集合Ω'の分散s2(j)(j=1,…,M)を各jについて求め、前記第j母集合Ωの分散σ2(j)(j=1,…,M)を用いて、
(n−1)×s2(j)/σ2(j)(j=1,…,M)を新たな第2確率変数χ2(j)(j=1,…,M)とする確率分布を生成し、該第2確率変数χ2(j)についての確率密度関数を用いて該第2確率変数χ2(j)(j=1,…,M)分布の裾の第2累積面積Pχ2(j)(j=1,…,M)を求めることから成ることとした。
上記モンテカルロ評価における危険不確定パラメータの検出方法では、上記構成とすることにより、上記確率Pz(j)と同様に、確率Pχ2(j)の値が小さい程、現実に起こった第j部分集合Ω'の分散s2(j)は母集団からn個を無作為に抽出した場合には起こりにくいことを示し、それに対応する不確定パラメータε(j)は失敗ケースに関わっている可能性が高いことになる。
(n−1)×s2(j)/σ2(j)(j=1,…,M)を新たな第2確率変数χ2(j)(j=1,…,M)とする確率分布を生成し、該第2確率変数χ2(j)についての確率密度関数を用いて該第2確率変数χ2(j)(j=1,…,M)分布の裾の第2累積面積Pχ2(j)(j=1,…,M)を求めることから成ることとした。
上記モンテカルロ評価における危険不確定パラメータの検出方法では、上記構成とすることにより、上記確率Pz(j)と同様に、確率Pχ2(j)の値が小さい程、現実に起こった第j部分集合Ω'の分散s2(j)は母集団からn個を無作為に抽出した場合には起こりにくいことを示し、それに対応する不確定パラメータε(j)は失敗ケースに関わっている可能性が高いことになる。
請求項4に記載のモンテカルロ評価における危険不確定パラメータの検出方法では、前記第1累積面積Pz(j)と前記第2累積面積Pχ2(j)とを各jについて比較し、小さい方を代表値P(j)(j=1,…,M)とし、該代表値P(j)を小さい順に並べることから成ることとした。
上記モンテカルロ評価における危険不確定パラメータの検出方法では、上記構成とすることにより、第1累積面積Pz(j)または第2累積面積Pχ2(j)に対応する不確定パラメータε(j)は失敗ケースに関わっている可能性が高い順番に並ぶことになる。
上記モンテカルロ評価における危険不確定パラメータの検出方法では、上記構成とすることにより、第1累積面積Pz(j)または第2累積面積Pχ2(j)に対応する不確定パラメータε(j)は失敗ケースに関わっている可能性が高い順番に並ぶことになる。
請求項5に記載のモンテカルロ評価における危険不確定パラメータの検出方法では、前記Mを用いて危険率αを1/(10×M)とし、前記代表値P(j)(j=1,…,M)の中の最小値と前記危険率αの大小を比較し、前記代表値P(j)が該危険率αより小さい場合はシステムの評価結果に重大な影響を及ぼす危険不確定パラメータと判定することから成ることとした。
上記モンテカルロ評価における危険不確定パラメータの検出方法では、不確定パラメータの数がM個であるとき、それぞれに対応する代表値P(j)もM個算出される。従って、M個算出されたP(j)の中の1個は偶然に1/M程度の確率であることは十分に起こり得る。このため、通常なら起こりえない程度の確率として、1/Mより一桁小さい値とした危険率αを1/(10×M)とすることにより好適に不確定パラメータを検定し、システムの評価結果に重大な影響を及ぼす危険不確定パラメータを検出することが出来る。
上記モンテカルロ評価における危険不確定パラメータの検出方法では、不確定パラメータの数がM個であるとき、それぞれに対応する代表値P(j)もM個算出される。従って、M個算出されたP(j)の中の1個は偶然に1/M程度の確率であることは十分に起こり得る。このため、通常なら起こりえない程度の確率として、1/Mより一桁小さい値とした危険率αを1/(10×M)とすることにより好適に不確定パラメータを検定し、システムの評価結果に重大な影響を及ぼす危険不確定パラメータを検出することが出来る。
本発明のモンテカルロ評価における危険不確定パラメータの検出方法によれば、失敗ケースに対応する不確定パラメータの確率分布に対し、全ケースに対応する不確定パラメータの確率分布に基づいて統計的検定を適用することにより、再度のモンテカルロ・シミュレーションを実施することなく、失敗ケースに関わった危険不確定パラメータを短時間で直ちに求めることが出来る。その結果、再度のモンテカルロ・シミュレーションに係る時間および作業が削減され、時間的コストが大幅に低減し且つ作業効率が大幅に向上するようになる。特に複数のパラメータの影響による失敗ケースが発生した場合であっても、再度のモンテカルロ・シミュレーションを実施することなく、失敗ケースに関わった危険不確定パラメータを短時間で直ちに求めることが出来る。
以下、図に示す実施の形態により本発明をさらに詳細に説明する。なお、これにより本発明が限定されるものではない。
図1は、本発明に係るモンテカルロ評価における危険不確定パラメータを検出する方法を示すフロー図である。
先ず、ステップS1では、失敗ケースに対応するn個の不確定パラメータ・ベクトル、ΕF 1,・・・,ΕF nを抽出する。ここで、不確定パラメータε(j)とは、センサ計測誤差、重心位置、温度、風、初期状態、システム特性など、その値がシステムの動作に影響を与えるが、その値が未知であるか又はその時々で値が異なるパラメータのことであり、不確定パラメータ・ベクトルとはこれらのパラメータを成分に持つベクトルのことである。
先ず、ステップS1では、失敗ケースに対応するn個の不確定パラメータ・ベクトル、ΕF 1,・・・,ΕF nを抽出する。ここで、不確定パラメータε(j)とは、センサ計測誤差、重心位置、温度、風、初期状態、システム特性など、その値がシステムの動作に影響を与えるが、その値が未知であるか又はその時々で値が異なるパラメータのことであり、不確定パラメータ・ベクトルとはこれらのパラメータを成分に持つベクトルのことである。
また、不確定パラメータ・ベクトルΕk(k=1,…,N)の構造を図2に示す。ここで、各成分εk(j)の右下の添え字kはシミュレーションの回数を示し、同( )内の数字jは、不確定パラメータの番号を示している。例えば、ε1(1)は1回目のシミュレーション時に乱数によって発生させたセンサ計測誤差であり、ε2(2)は2回目のシミュレーション時に乱数によって発生させた重心位置であることを示す。
ステップS2では、図3に示すように、第j部分集合Ω'の確率変数εF k(j)(k=1,…,n)の平均ε'F(j)および分散s2(j)を求める。なお、平均および分散は下記式により求める。
ε'F(j)=(1/n)×Σi=1 nεF i(j)(j=1,…,M)
s2(j)=1/(n−1)×Σi=1 n{εF i(j)−ε'F(j)}2(j=1,…,M)
ε'F(j)=(1/n)×Σi=1 nεF i(j)(j=1,…,M)
s2(j)=1/(n−1)×Σi=1 n{εF i(j)−ε'F(j)}2(j=1,…,M)
次に、ステップS3では、第j母集合Ωの確率変数εk(j)(k=1,…,N)の平均μ(j)および分散σ2(j)を用いて、第j部分集合Ω'の平均ε'F(j)を下記の通り標準化し、新たな確率変数z(j)(j=1,…,M)を生成する。
z(j)=(ε’F(j)−μ(j))/(σ2(j)/n)1/2(j=1,…,M)
また、第j部分集合Ω'の分散s2(j)についても、第j母集合Ωの確率変数εk(j)(k=1,…,N)の平均μ(j)および分散σ2(j)を用いて、新たな確率変数χ2(j)(j=1,…,M)を生成する。
χ2(j)=(n−1)×s2(j)/σ2(j)(j=1,…,M)
このとき、z(j)はN(0,1)の標準正規分布に従い、χ2(j)は自由度(n−1)のχ2分布に従う。
z(j)=(ε’F(j)−μ(j))/(σ2(j)/n)1/2(j=1,…,M)
また、第j部分集合Ω'の分散s2(j)についても、第j母集合Ωの確率変数εk(j)(k=1,…,N)の平均μ(j)および分散σ2(j)を用いて、新たな確率変数χ2(j)(j=1,…,M)を生成する。
χ2(j)=(n−1)×s2(j)/σ2(j)(j=1,…,M)
このとき、z(j)はN(0,1)の標準正規分布に従い、χ2(j)は自由度(n−1)のχ2分布に従う。
なお、図4は、本発明に係る不確定パラメータ・ベクトルの検出方法の中核を成す抽出原理を示す説明図である。
例えば、全ケースに係る第j母集合Ωの確率分布と、失敗ケースに係る第j部分集合の確率分布とを比較した場合、母集合から無作為に抽出したとは考えられない程、第j部分集合の分布が異なっている。一般に、母集団から無作為に抽出された標本の分布は、母集団の分布に近似するようになる。従って、この不確定パラメータ、例えば第j番目の不確定パラメータε(j)は失敗ケースに関わっている可能性が十分考えられる。そして、どの程度異なっていれば、その不確定パラメータは失敗ケースに関わっていると判定してよいかという具体的手順については、後述する各不確定パラメータに対応するM個の累積確率P(j)を危険率αで検定することにより行う。
例えば、全ケースに係る第j母集合Ωの確率分布と、失敗ケースに係る第j部分集合の確率分布とを比較した場合、母集合から無作為に抽出したとは考えられない程、第j部分集合の分布が異なっている。一般に、母集団から無作為に抽出された標本の分布は、母集団の分布に近似するようになる。従って、この不確定パラメータ、例えば第j番目の不確定パラメータε(j)は失敗ケースに関わっている可能性が十分考えられる。そして、どの程度異なっていれば、その不確定パラメータは失敗ケースに関わっていると判定してよいかという具体的手順については、後述する各不確定パラメータに対応するM個の累積確率P(j)を危険率αで検定することにより行う。
次いで、ステップS4では、z(j)分布の裾の累積面積(累積確率)Pz(j)を求める。ここで、標準正規分布の確率密度関数を図5に示す。z(j)分布の裾の面積を求めるとは図中の斜線部分面線を求めることであるが、z(j)の分布は0を中心として左右対称であるため、z≧|z(j)|における面積(確率)をPz(j)とする。Pz(j)の値が小さい程、現実に起こった平均ε'F(j)は、第j母集合Ωからn個を無作為に抽出した場合には起こりにくいことを示している。つまり、累積確率Pz(j)の値が小さい程、無作為に抽出したとは考えにくく、失敗ケースに関わっている可能性が高いことになる。
次いで、χ2(j)分布の裾の累積面積(累積確率)Pχ2(j)を求める。ここで、図6にχ2分布を示す。χ2の値は常に0以上であり、左右対称でもない。従って、先ず下側累積確率Pχ2(j)を求める。Pχ2(j)≦0.5ならばPχ2(j)が分布の裾の面積となるが、Pχ2(j)>0.5のときは上側累積確率が裾の面積となるため、Pχ2(j)←1−Pχ2(j)とする。このように定義することにより、Pz(j)と同様に、Pχ2(j)の値が小さい程、失敗ケースに関わっている可能性が高くなる。
次いで、ステップS5では、各不確定パラメータε(j)についての累積確率値P(j)(j=1,…,M)の値を、下記に示すように、z(j)の累積確率Pz(j)と、χ2(j)の累積確率Pχ2(j)とを比較して、小さい方の値として定義する。
P(j)=min{Pz(j),Pχ2(j)}
P(j)=min{Pz(j),Pχ2(j)}
次いで、ステップS6では、累積確率値P(j)(j=1,…,M)を小さい順に並べ替える。この順番に対応する不確定パラメータε(j)(j=1,…,n)は失敗ケースに関わっている可能性が高い順番に並ぶ。
次いで、ステップS7では、危険率の目安を設定する。統計的検定においては通常危険率αを設定し、分布の裾の面積(累積確率)がαよりも小さい時、帰無仮説を棄却して対応する不確定パラメータは失敗ケースに影響していると判断する。この場合、累積確率値P(j)の小さい順に失敗ケースに関わっている可能性が高いが、累積確率値P(j)の値がどの程度なら対応する不確定パラメータε(j)(j=1,…,M)が失敗ケースに関わっていると判断して良いか否かの目安(危険率)を設定する。
不確定パラメータε(j)の数がM個であるとき、それぞれに対応する累積確率P(j)もM個算出される。したがって、M個算出されたP(j)の中の1個は偶然に1/M程度の確率であることは十分に起こり得る。このため、通常なら起こりえない程度の確率として、1/Mよりも一桁小さい値として次式により危険率を設定する。
α=1/(10×M)
不確定パラメータε(j)の数がM個であるとき、それぞれに対応する累積確率P(j)もM個算出される。したがって、M個算出されたP(j)の中の1個は偶然に1/M程度の確率であることは十分に起こり得る。このため、通常なら起こりえない程度の確率として、1/Mよりも一桁小さい値として次式により危険率を設定する。
α=1/(10×M)
最後に、ステップS8では、ステップS6で得られたM個の累積確率P(j)を危険率αで検定する。すなわち、累積確率P(j)がこの値よりも小さいときは、対応する不確定パラメータε(j)は失敗ケースに影響している可能性が高いと判定する。
以上のことをまとめると、図7から図9のようになる。
以上のことをまとめると、図7から図9のようになる。
本発明のモンテカルロ評価における危険不確定パラメータの検出方法を、自動着陸実験のモンテカルロ評価結果に適用した結果を示す。なお、自動着陸実験(ALFLEX)は、宇宙往還機の自動着陸技術を確立するために1996年に実施されたものである。実験機はヘリコプタに吊されて飛行し、その後滑走路手前2700m、高度1500mにおいて分離される。その後実験機は各種センサ信号を搭載計算機に取り込み、誘導制御ロジックに従ってコントロールされ、自動飛行により滑走路に着陸する。また、実験の概要を図10に示す。
[モンテカルロ評価]
本実験では自動着陸実験前に試験飛行による動作確認を実施することは困難であり、初飛行時に自動着陸を成功させなければならない。従って、システムの事前評価が非常に重要となる。事前評価のためにシステムの数学モデルを構築するが、現実にはセンサ計測誤差、空気力の誤差、質量や慣性モーメントの誤差、風条件、初期条件等の様々な不確定要因(不確定パラメータ)が存在し、数学モデルでは記述できない要素が多く存在する。モンテカルロ・シミュレーションは、これらの影響も含めた評価が可能な手段であり、本実験の事前評価手法として最適である。
モンテカルロ評価においては、各不確定パラメータε(j)の分布を仮定し、その分布に基づき乱数を用いて不確定パラメータεk(j)(j=1,…,M、k=1,…,N)の値を発生させる。これらを同時に加え、多数回のシミュレーションを実行することによってミッション達成確率が推定できる。本実験のモンテカルロ評価において、加えられた不確定パラメータの数は100(M=100)であり、その項目を図11に示す。
本実験では自動着陸実験前に試験飛行による動作確認を実施することは困難であり、初飛行時に自動着陸を成功させなければならない。従って、システムの事前評価が非常に重要となる。事前評価のためにシステムの数学モデルを構築するが、現実にはセンサ計測誤差、空気力の誤差、質量や慣性モーメントの誤差、風条件、初期条件等の様々な不確定要因(不確定パラメータ)が存在し、数学モデルでは記述できない要素が多く存在する。モンテカルロ・シミュレーションは、これらの影響も含めた評価が可能な手段であり、本実験の事前評価手法として最適である。
モンテカルロ評価においては、各不確定パラメータε(j)の分布を仮定し、その分布に基づき乱数を用いて不確定パラメータεk(j)(j=1,…,M、k=1,…,N)の値を発生させる。これらを同時に加え、多数回のシミュレーションを実行することによってミッション達成確率が推定できる。本実験のモンテカルロ評価において、加えられた不確定パラメータの数は100(M=100)であり、その項目を図11に示す。
[評価結果]
100個の不確定パラメータε(j)(j=1,…,100)を同時に加えた飛行シミュレーションを、それらの値を更新しながら5000回(N=5000)実行した。その結果、途中で姿勢運動が不安定となり、飛行不能となって着陸まで至らなかったケースが17回(n=17)現れた。各不確定パラメータε(j)を一つずつ加えた飛行シミュレーション(感度解析)では、飛行不能となるケースは発生していない。このため、17回の飛行不能となったケースは複数の不確定パラメータε(j)が複合して影響しているものと考えられる。そこで、次に、本発明のモンテカルロ評価における不確定パラメータの検出方法を用いて、着陸失敗の原因となる危険不確定パラメータε(j)を特定することを試みる。
100個の不確定パラメータε(j)(j=1,…,100)を同時に加えた飛行シミュレーションを、それらの値を更新しながら5000回(N=5000)実行した。その結果、途中で姿勢運動が不安定となり、飛行不能となって着陸まで至らなかったケースが17回(n=17)現れた。各不確定パラメータε(j)を一つずつ加えた飛行シミュレーション(感度解析)では、飛行不能となるケースは発生していない。このため、17回の飛行不能となったケースは複数の不確定パラメータε(j)が複合して影響しているものと考えられる。そこで、次に、本発明のモンテカルロ評価における不確定パラメータの検出方法を用いて、着陸失敗の原因となる危険不確定パラメータε(j)を特定することを試みる。
[不確定パラメータの検出]
着陸不能となった17回の飛行ケースに対応する17個の不確定パラメータ・ベクトルΕF 1,…,ΕF 17を取り出し、これらの不確定パラメータ・ベクトルに対し本発明のモンテカルロ評価における危険不確定パラメータの検出方法を適用した。100個の不確定パラメータε(j)のそれぞれについて、対応する第j部分集合Ω'の平均ε'F(j)および分散s2(j)を求め、これらの平均および分散に対して、上記ステップS3にて示した処理を行って新たな確率変数z(j)およびχ2(j)を生成し、そのz(j)分布から対応する裾部の累積確率Pz(j)を,そのχ2(j)分布から対応する裾部の累積確率Pχ2(j)を各々求め、上記ステップS5にて示した処理を行って累積確率P(j)(j=1,…,100)を算出し、このP(j)(j=1,…,100)の小さい順に対応する不確定パラメータを並べた。また、このモンテカルロ評価で加えた不確定パラメータは100個であるので、着陸不能ケースに影響があると判定するP(j)値の目安(危険率)として、
α=1/(10・100)=1/1000=1×10-3
とした。従って、P(j)値がα以下であるときは、対応する不確定パラメータは着陸不能ケースに影響していると考えられる。α以下となったP(j)値と対応する不確定パラメータを表1に示す。
表1から、着陸不能ケースに影響していると考えられる危険不確定パラメータは、エレボン・レートリミッタ、ラダー・レートリミッタ、分離時のヨーレートであり、これらに対応するP(j)値は確かにαよりもかなり小さい。図12に、これらの検出された3つの不確定パラメータに関する17個のデータの度数分布を示す。
エレボン・レートリミッタ及びラダー・レートリミッタは、舵のアクチュエータの回転角速度を表す。また、表1から、これらについてはPz(j)の値がP(j)値となっており、着陸に失敗した17回の平均ε'F(j)はノミナル状態に比べて小さくなっていることが図12からもわかる。つまり、ノミナル状態に比べて舵の回転角速度がかなり小さくなっていることを意味し、制御指令に舵が追従できていない可能性があることを示す。また、分離時ヨーレートは、実験機をヘリコプタから分離する時に、ヨー方向の角速度が発生していることを表す。
分離時ヨーレートについては、Pχ2(j)の値がP(j)値となっており、着陸失敗となった17回の分散がノミナル状態よりも大きいことが図12からもわかる。これは正負に依らずヨーレートの絶対値が大きい場合を表し、分離時にヨー方向の運動が励起されている状態である。
着陸不能となった17回の飛行ケースに対応する17個の不確定パラメータ・ベクトルΕF 1,…,ΕF 17を取り出し、これらの不確定パラメータ・ベクトルに対し本発明のモンテカルロ評価における危険不確定パラメータの検出方法を適用した。100個の不確定パラメータε(j)のそれぞれについて、対応する第j部分集合Ω'の平均ε'F(j)および分散s2(j)を求め、これらの平均および分散に対して、上記ステップS3にて示した処理を行って新たな確率変数z(j)およびχ2(j)を生成し、そのz(j)分布から対応する裾部の累積確率Pz(j)を,そのχ2(j)分布から対応する裾部の累積確率Pχ2(j)を各々求め、上記ステップS5にて示した処理を行って累積確率P(j)(j=1,…,100)を算出し、このP(j)(j=1,…,100)の小さい順に対応する不確定パラメータを並べた。また、このモンテカルロ評価で加えた不確定パラメータは100個であるので、着陸不能ケースに影響があると判定するP(j)値の目安(危険率)として、
α=1/(10・100)=1/1000=1×10-3
とした。従って、P(j)値がα以下であるときは、対応する不確定パラメータは着陸不能ケースに影響していると考えられる。α以下となったP(j)値と対応する不確定パラメータを表1に示す。
表1から、着陸不能ケースに影響していると考えられる危険不確定パラメータは、エレボン・レートリミッタ、ラダー・レートリミッタ、分離時のヨーレートであり、これらに対応するP(j)値は確かにαよりもかなり小さい。図12に、これらの検出された3つの不確定パラメータに関する17個のデータの度数分布を示す。
エレボン・レートリミッタ及びラダー・レートリミッタは、舵のアクチュエータの回転角速度を表す。また、表1から、これらについてはPz(j)の値がP(j)値となっており、着陸に失敗した17回の平均ε'F(j)はノミナル状態に比べて小さくなっていることが図12からもわかる。つまり、ノミナル状態に比べて舵の回転角速度がかなり小さくなっていることを意味し、制御指令に舵が追従できていない可能性があることを示す。また、分離時ヨーレートは、実験機をヘリコプタから分離する時に、ヨー方向の角速度が発生していることを表す。
分離時ヨーレートについては、Pχ2(j)の値がP(j)値となっており、着陸失敗となった17回の分散がノミナル状態よりも大きいことが図12からもわかる。これは正負に依らずヨーレートの絶対値が大きい場合を表し、分離時にヨー方向の運動が励起されている状態である。
[検出結果の確認]
次に、これらの検出された3つの不確定パラメータ(エレボン・レートリミッタ,ラダー・レートリミッタ,分離時のヨーレート)が、本当に17回の着陸失敗ケースの原因であるのかどうかを確認した。このことを確認するために、表1に示される3つの不確定パラメータを除き、他の条件は全て同一として、再度5000回のモンテカルロ評価を実施した。その結果、着陸失敗となるケースは一回も発生しなかった。
また念のため、検出された3つの不確定パラメータを一つずつ除いた場合、2つだけ除いた場合について、それぞれ5000回のモンテカルロ評価を実施したが、いずれの場合も着陸失敗は10回以下となった。以上のことから、検出された3つの不確定パラメータはいずれも着陸失敗に影響していることが確認できた。
解析の結果、実験機の分離時にヨーレートが発生しており、かつ、エレボン・レートリミッタ又はラダー・レートリミッタが小さくなっているときには、機体のヨー運動を抑えるための制御指令に舵を動かすアクチュエータが追従できていないことがわかった。結果として、姿勢制御不能となり着陸する前に姿勢運動が発散する現象が発生している。なお、図13に参考として、ノミナル状態の機体モデルに「エレボン・レートリミッタ+分離時ヨーレート」、「ラダー・レートリミッタ+分離時ヨーレート」の不確定パラメータの組み合わせで飛行シミュレーションを実施した結果を示す。両ケース共に、実験機分離後直ちに姿勢運動が不安定となっており、10秒も経たないうちに異常な姿勢となって飛行不能となっていることがわかる。
以上より、本発明によるモンテカルロ評価における危険不確定パラメータの検出方法によれば、不確定パラメータ単独ではシステムに対し影響を及ぼさなくても、複数のパラメータが複合的に作用し合うことによりシステムに対し影響を及ぼす場合であっても、失敗の原因となる危険不確定パラメータは確実に検出される。
次に、これらの検出された3つの不確定パラメータ(エレボン・レートリミッタ,ラダー・レートリミッタ,分離時のヨーレート)が、本当に17回の着陸失敗ケースの原因であるのかどうかを確認した。このことを確認するために、表1に示される3つの不確定パラメータを除き、他の条件は全て同一として、再度5000回のモンテカルロ評価を実施した。その結果、着陸失敗となるケースは一回も発生しなかった。
また念のため、検出された3つの不確定パラメータを一つずつ除いた場合、2つだけ除いた場合について、それぞれ5000回のモンテカルロ評価を実施したが、いずれの場合も着陸失敗は10回以下となった。以上のことから、検出された3つの不確定パラメータはいずれも着陸失敗に影響していることが確認できた。
解析の結果、実験機の分離時にヨーレートが発生しており、かつ、エレボン・レートリミッタ又はラダー・レートリミッタが小さくなっているときには、機体のヨー運動を抑えるための制御指令に舵を動かすアクチュエータが追従できていないことがわかった。結果として、姿勢制御不能となり着陸する前に姿勢運動が発散する現象が発生している。なお、図13に参考として、ノミナル状態の機体モデルに「エレボン・レートリミッタ+分離時ヨーレート」、「ラダー・レートリミッタ+分離時ヨーレート」の不確定パラメータの組み合わせで飛行シミュレーションを実施した結果を示す。両ケース共に、実験機分離後直ちに姿勢運動が不安定となっており、10秒も経たないうちに異常な姿勢となって飛行不能となっていることがわかる。
以上より、本発明によるモンテカルロ評価における危険不確定パラメータの検出方法によれば、不確定パラメータ単独ではシステムに対し影響を及ぼさなくても、複数のパラメータが複合的に作用し合うことによりシステムに対し影響を及ぼす場合であっても、失敗の原因となる危険不確定パラメータは確実に検出される。
なお、上記実施例では、本発明をモンテカルロ法に係る評価結果に対し適用したが、モンテカルロ法に限らず、ある確率分布に従う複数の不確定パラメータを設定してシステムの評価を行うシミュレーション技法におけるシステムの評価に重大な影響を及ぼす危険不確定パラメータの検出に対し好適に適用することが出来る。
本発明のモンテカルロ評価における危険不確定パラメータの検出方法は、実環境での試験が容易ではなく、かつ、環境条件・システムの特性・初期条件等様々な不確定パラメータの存在下で動作することが要求されるシステムに対して、事前のモンテカルロ評価が有効であるシステム、例えば航空宇宙機・船舶・各種プラント等のシステムの評価に対し利用可能である。
Εk 不確定パラメータ・ベクトル
ΕF k 失敗例に対応する不確定パラメータ・ベクトル
μ(j) 第j母集合における不確定パラメータの平均
σ2(j) 第j母集合における不確定パラメータの分散
ε'F(j) 第j部分集合における不確定パラメータの平均
s2(j) 第j部分集合における不確定パラメータの分散
z(j) ε'F(j)を標準化した確率変数
χ2(j) s2(j)を標準化した確率変数
Pz(j) z(j)分布の裾の累積面積
Pχ2(j) χ2(j)の裾の累積面積
P(j) 裾の累積面積の代表値
100 モンテカルロ評価における危険不確定パラメータの検出フロー
ΕF k 失敗例に対応する不確定パラメータ・ベクトル
μ(j) 第j母集合における不確定パラメータの平均
σ2(j) 第j母集合における不確定パラメータの分散
ε'F(j) 第j部分集合における不確定パラメータの平均
s2(j) 第j部分集合における不確定パラメータの分散
z(j) ε'F(j)を標準化した確率変数
χ2(j) s2(j)を標準化した確率変数
Pz(j) z(j)分布の裾の累積面積
Pχ2(j) χ2(j)の裾の累積面積
P(j) 裾の累積面積の代表値
100 モンテカルロ評価における危険不確定パラメータの検出フロー
Claims (5)
- モンテカルロ法を用いたシミュレーションにおいてその評価結果に重大な影響を及ぼす危険不確定パラメータを検出する方法であって、乱数によって発生させるM個の不確定パラメータの内で第j不確定パラメータをε(j)(j=1,…,M)と、全シミュレーションの回数をNと、第k回目のシミュレーションにおいて発生させたM個の不確定パラメータの内で第j不確定パラメータをεk(j)(j=1,…,M)とする時、該不確定パラメータを成分に持つ不確定パラメータ・ベクトルΕk=t[εk(j)|j=1,…,M](k=1,…,N)を要素に持つ母集合[Εk|k=1,…,N]から、失敗例に対応するn個の不確定パラメータ・ベクトルΕF k(k=1,…,n)を抽出し、該抽出した各不確定パラメータ・ベクトルの第j成分εF k(j)(k=1,…,n)のみから構成される第j部分集合Ω'=[εF k(j)|k=1,…,n]の平均ε'F(j)および分散s2(j)を新たな確率変数とする確率分布を求め、前記母集合の第j不確定パラメータεk(j)(k=1,…,N)のみから構成される第j母集合Ω=[εk(j)|k=1,…,N]の確率分布に基づいて、前記第j部分集合Ω'に係る確率変数の分布に対し統計的検定を適用することにより評価結果に重大な影響を及ぼす危険不確定パラメータを求めることを特徴とするモンテカルロ評価における危険不確定パラメータの検出方法。
- 前記統計的検定は、前記第j部分集合Ω'の平均ε'F(j)(j=1,…,M)を各jについて求め、次いで前記第j母集合Ωの平均μ(j)(j=1,…,M)および分散σ2(j)(j=1,…,M)を用いて、
(ε'F(j)−μ(j))/(σ2(j)/n)1/2(j=1,…,M)を新たな第1確率変数z(j)(j=1,…,M)とする確率分布を生成し、該第1確率変数z(j)についての確率密度関数を用いて該第1確率変数z(j)分布の裾の第1累積面積Pz(j)(j=1,…,M)を求めることから成る請求項1に記載のモンテカルロ評価における危険不確定パラメータの検出方法。 - 前記統計的検定は、前記第j部分集合Ω'の分散s2(j)(j=1,…,M)を各jについて求め、前記第j母集合Ωの分散σ2(j)(j=1,…,M)を用いて、
(n−1)×s2(j)/σ2(j)(j=1,…,M)を新たな第2確率変数χ2(j)(j=1,…,M)とする確率分布を生成し、該第2確率変数χ2(j)についての確率密度関数を用いて該第2確率変数χ2(j)(j=1,…,M)分布の裾の第2累積面積Pχ2(j)(j=1,…,M)を求めることから成る請求項1に記載のモンテカルロ評価における危険不確定パラメータの検出方法。 - 前記第1累積面積Pz(j)と前記第2累積面積Pχ2(j)とを各jについて比較し、小さい方を代表値P(j)(j=1,…,M)とし、該代表値P(j)を小さい順に並べることから成る請求項1から3の何れかに記載のモンテカルロ評価における危険不確定パラメータの検出方法。
- 前記Mを用いて危険率αを1/(10×M)とし、前記代表値P(j)(j=1,…,M)の中の最小値と前記危険率αの大小を比較し、前記代表値P(j)が該危険率αより小さい場合はシステムの評価結果に重大な影響を及ぼす危険不確定パラメータと判定することから成る請求項4に記載のモンテカルロ評価における危険不確定パラメータの検出方法。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2006030256A JP2007213172A (ja) | 2006-02-07 | 2006-02-07 | モンテカルロ評価における危険不確定パラメータの検出方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2006030256A JP2007213172A (ja) | 2006-02-07 | 2006-02-07 | モンテカルロ評価における危険不確定パラメータの検出方法 |
Publications (1)
Publication Number | Publication Date |
---|---|
JP2007213172A true JP2007213172A (ja) | 2007-08-23 |
Family
ID=38491569
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP2006030256A Pending JP2007213172A (ja) | 2006-02-07 | 2006-02-07 | モンテカルロ評価における危険不確定パラメータの検出方法 |
Country Status (1)
Country | Link |
---|---|
JP (1) | JP2007213172A (ja) |
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2016030938A1 (ja) * | 2014-08-25 | 2016-03-03 | 株式会社日立製作所 | シミュレーションシステム、及びシミュレーション方法 |
CN109871515A (zh) * | 2018-11-22 | 2019-06-11 | 常州天正工业发展股份有限公司 | 一种根据设备数据评估危险发生概率的系统和算法 |
CN111125876A (zh) * | 2019-11-18 | 2020-05-08 | 中国空间技术研究院 | 一种无铅元器件互连焊点锡须触碰风险评估方法 |
CN114329974A (zh) * | 2021-12-29 | 2022-04-12 | 武汉大学 | 基于蒙特卡洛模拟的城市供水管网地震损伤评估方法 |
-
2006
- 2006-02-07 JP JP2006030256A patent/JP2007213172A/ja active Pending
Cited By (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2016030938A1 (ja) * | 2014-08-25 | 2016-03-03 | 株式会社日立製作所 | シミュレーションシステム、及びシミュレーション方法 |
JPWO2016030938A1 (ja) * | 2014-08-25 | 2017-06-15 | 株式会社日立製作所 | シミュレーションシステム、及びシミュレーション方法 |
CN109871515A (zh) * | 2018-11-22 | 2019-06-11 | 常州天正工业发展股份有限公司 | 一种根据设备数据评估危险发生概率的系统和算法 |
CN109871515B (zh) * | 2018-11-22 | 2023-07-07 | 常州天正工业发展股份有限公司 | 一种根据设备数据评估危险发生概率的系统和算法 |
CN111125876A (zh) * | 2019-11-18 | 2020-05-08 | 中国空间技术研究院 | 一种无铅元器件互连焊点锡须触碰风险评估方法 |
CN111125876B (zh) * | 2019-11-18 | 2024-02-09 | 中国空间技术研究院 | 一种无铅元器件互连焊点锡须触碰风险评估方法 |
CN114329974A (zh) * | 2021-12-29 | 2022-04-12 | 武汉大学 | 基于蒙特卡洛模拟的城市供水管网地震损伤评估方法 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Abbaspour et al. | Neural adaptive observer-based sensor and actuator fault detection in nonlinear systems: Application in UAV | |
Samy et al. | Survey and application of sensor fault detection and isolation schemes | |
Odendaal et al. | Actuator fault detection and isolation: An optimised parity space approach | |
Fravolini et al. | Data-driven schemes for robust fault detection of air data system sensors | |
Aboutalebi et al. | A novel sensor fault detection in an unmanned quadrotor based on adaptive neural observer | |
Zhang et al. | Aircraft damage identification and classification for database-driven online flight-envelope prediction | |
Rotondo et al. | Diagnosis of icing and actuator faults in UAVs using LPV unknown input observers | |
Zolghadri | The challenge of advanced model-based FDIR for real-world flight-critical applications | |
Cimpoesu et al. | Fault detection and diagnosis using parameter estimation with recursive least squares | |
Xu et al. | Fault tolerant control for near space vehicle: a survey and some new results | |
Samy et al. | Sensor fault detection and accommodation using neural networks with application to a non-linear unmanned air vehicle model | |
JP2007213172A (ja) | モンテカルロ評価における危険不確定パラメータの検出方法 | |
Li et al. | Operational anomaly detection in flight data using a multivariate gaussian mixture model | |
Al-Haddad et al. | An intelligent quadcopter unbalance classification method based on stochastic gradient descent logistic regression | |
Perhinschi et al. | Development of immunity-based framework for aircraft abnormal conditions detection, identification, evaluation, and accommodation | |
Azam et al. | In-flight fault detection and isolation in aircraft flight control systems | |
Gururajan et al. | Evaluation of sensor failure detection, identification and accommodation (SFDIA) performance following common-mode failures of pitot tubes | |
Venkataraman et al. | Model-based detection and isolation of rudder faults for a small UAS | |
Samy et al. | Fault detection and flight data measurement: Demonstrated on unmanned air vehicles using neural networks | |
Hansen et al. | Control surface fault diagnosis with specified detection probability—Real event experiences | |
Jing et al. | Fault detection in Quadrotor MAV | |
Moncayo et al. | Aircraft failure detection and identification over an extended flight envelope using an artificial immune system | |
Lu | Fault diagnosis and fault-tolerant control for aircraft subjected to sensor and actuator faults | |
Prabhakaran et al. | System safety analysis for critical system applications using Bayesian networks | |
Guven et al. | Two-Stage Kalman Filter Based Estimation of Boeing 747 Actuator/Control Surface Stuck Faults |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20080702 |
|
A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20080829 |
|
A02 | Decision of refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A02 Effective date: 20080924 |