(57)【特許請求の範囲】
【請求項１】
ベクトル写像fを用いてm (m≧1)次元の非有界領域M = (-∞,+∞)×…×(-∞,+∞)に対してm次元被積分関数Aの多重積分を計算する計算装置であって、
当該ベクトル写像fは、実数の成分からなるm次元のベクトルを実数の成分からなるm次元のベクトルに変換し、
当該ベクトル写像fを任意のｍ次元のベクトルx∈Mに適用した結果f(x)は当該ベクトルxから計算可能であり、
当該ベクトル写像fを所定のm次元のベクトルuに繰り返し適用した結果u，f(u)，f(f(u))，f(f(f(u)))，… の極限分布の多次元密度関数ρを任意のｍ次元のベクトルx∈Mに適用した結果ρ(x)は当該ベクトルxから計算可能であり、
当該計算装置は、
m次元のベクトルx=(x1,x2,…,xm)を記憶する第１の記憶部と、
スカラー値wを記憶する第２の記憶部と、
前記第１の記憶部に記憶されたベクトルxに、前記ベクトル写像fを適用した結果のベクトルx'=f(x)=(x'1,x'2,…,x'm)を計算する第１の計算部と、
前記第１の記憶部に記憶されたベクトルxから、スカラー値w'=A(x)/ρ(x)を計算する第２の計算部と、
前記第１の計算部により計算されたベクトルx'を前記第１の記憶部に記憶されるベクトルxに代入して、前記第１の記憶部に記憶されるベクトルxを更新し、
前記第２の計算部により計算されたスカラー値w'を前記第２の記憶部に記憶されるスカラー値wに加算して、前記第２の記憶部に記憶されるスカラー値wを更新する更新部と、
前記更新部による更新回数がc (c≧1)回となった場合、前記第２の記憶部に記憶されたスカラー値wから、スカラー値s=w/(c+1)を計算して、当該スカラー値sを多重積分の結果として出力する出力部と、
を備え、
前記第２の記憶部に最初に記憶されるスカラー値は、前記第１の記憶部に最初に記憶されるm次元のベクトルに関数Aを適用したものを、前記第１の記憶部に最初に記憶されるm次元のベクトルに密度関数ρを適用したもので除算した結果である
ことを特徴とする計算装置。
【請求項２】
ベクトル写像fを用いて、m (m≧1)次元の非有界領域M = (-∞,+∞)×…×(-∞,+∞)に対して多次元被積分関数Aの多重積分を計算する計算装置であって、
当該ベクトル写像fは、実数の成分からなるm次元のベクトルを実数の成分からなるm次元のベクトルに変換し、
当該ベクトル写像fを任意のｍ次元のベクトルx∈Mに適用した結果f(x)は当該ベクトルxから計算可能であり、
当該ベクトル写像fを所定のm次元のベクトルuに繰り返し適用した結果u，f(u)，f(f(u))，f(f(f(u)))，… の極限分布の多次元密度関数ρを任意のｍ次元のベクトルx∈Mに適用した結果ρ(x)は当該ベクトルxから計算可能であり、
当該計算装置は、
m次元のベクトルx=(x1,x2,…,xm)を記憶する第１の記憶部と、
スカラー値wを記憶する第２の記憶部と、
前記第１の記憶部に記憶されたベクトルxに、前記ベクトル写像fを適用した結果のベクトルx'=f(x)=(x'1,x'2,…,x'm)を計算する第１の計算部と、
前記第１の記憶部に記憶されたベクトルxから、スカラー値w'=A(x)/ρ(x)を計算する第２の計算部と、
前記第１の計算部により計算されたベクトルx'を前記第１の記憶部に記憶されるベクトルxに代入して、前記第１の記憶部に記憶されるベクトルxを更新し、
前記第２の計算部により計算されたスカラー値w'を前記第２の記憶部に記憶されるスカラー値wに加算して、前記第２の記憶部に記憶されるスカラー値wを更新する更新部と、
前記更新部による更新回数がc (c≧1)回となった場合、前記第２の記憶部に記憶されたスカラー値wから、スカラー値s'=w/cを計算して、当該スカラー値s'を多重積分の結果として出力する出力部と、
を備え、
前記第２の記憶部に最初に記憶されるスカラー値は０である
ことを特徴とする計算装置。
【請求項３】
前記ベクトル写像fは、m次元ベクトルu = (u1，u2，…，um)とm個の写像f1，f2，…，fmとにより、
f(u) = (f1(u1)，f2(u2)，…，fm(um)))
と定義され、当該m個の写像f1，f2，…，fmはいずれもカオス写像である
ことを特徴とする請求項１または２に記載の計算装置。
【請求項４】
前記ベクトル写像fとして複数のベクトル写像g1，g2，…，gk (k≧2)を用意し、そのそれぞれについて、前記更新回数cを変化させて前記出力部が順次出力するスカラー値の収束速度を取得する収束速度取得部と、
前記収束速度取得部により取得された収束速度が最も速いベクトル写像gh (1≦h≦k)を選択する選択部と、
前記ベクトル写像fとして前記ベクトル写像ghを用いて、前記更新回数cにかえて更新回数c' (c'>c)を用いて前記出力部にスカラー値を出力させる出力制御部と、
をさらに備えることを特徴とする請求項１から３のいずれか１項に記載の計算装置。
【請求項５】
前記ベクトル写像fを所定のm次元のベクトルu=(u1,u2,…,um)に０回以上適用して得られるベクトル列
u，f(u)，f(f(u))，f(f(f(u)))，…
の極限分布の多次元密度関数ρは、性質
ρ(u) = Πi=0m ρi(ui)；
ρi(ui) 〜 c-i|ui|-(1+a)   for ui → -∞；
ρi(ui) 〜 c+i|ui|-(1+a)   for ui → +∞；
(a>0，1≦i≦m，c-i>0，c+i>0)
を満たす
ことを特徴とする請求項１から３のいずれか１項に記載の計算装置。
【請求項６】
前記ベクトル写像ｆは、1≦i≦mについて定義される関数
fi(t) = gi(di t)/di   (di>0)
により
f(u)=(f1(u1)，f2(u2)，…，fm(um))
と定義され、写像giは、以下の写像φj (1≦j≦8)のいずれかと自然数ni (ni≧2)とにより、
gi(φj(θ)) = φj(niθ)；
φ1(θ) = -sgn(tanθ)/|tanθ|1/a；
φ2(θ) = -sgn(tanθ)×|tanθ|1/a；
φ3(θ) = -sgn(cosθ)/|tanθ|1/a；
φ4(θ) = -sgn(cosθ)×|tanθ|1/a；
φ5(θ) =  sgn(cosθ)/|tanθ|1/a；
φ6(θ) =  sgn(cosθ)×|tanθ|1/a；
φ7(θ) =  sgn(sinθ)/|tanθ|1/a；
φ8(θ) =  sgn(sinθ)×|tanθ|1/a；
sgn(t) =  1   for t>0；
sgn(t) =  0   for t=0；
sgn(t) = -1   for t<0
と定義される
ことを特徴とする請求項５に記載の計算装置。
【請求項７】
前記定数aとして複数の正数q1，q2，…，qk (k≧2)を用意して、そのそれぞれについて前記写像fを定義し、前記更新回数cを変化させて前記出力部が順次出力するスカラー値の収束速度を取得する収束速度取得部と、
前記収束速度取得部により取得された収束速度が最も速い正数qh (1≦h≦k)を選択する選択部と、
前記定数aとして前記正数qhを用いて前記写像fを定義し、前記更新回数cにかえて更新回数c' (c'>c)を用いて前記出力部にスカラー値を出力させる出力制御部と、
をさらに備えることを特徴とする請求項５または６に記載の計算装置。
【請求項８】
前記写像giを前記写像φjのいずれか複数を用いて定義し、そのそれぞれについて、前記更新回数cを変化させて前記出力部が順次出力するスカラー値の収束速度を取得する収束速度取得部と、
前記収束速度取得部により取得された収束速度が最も速い写像φjを選択する選択部と、
前記写像giを前記選択部により選択された写像φjを用いて定義し、前記更新回数cにかえて更新回数c' (c'>c)を用いて前記出力部にスカラー値を出力させる出力制御部と、
をさらに備えることを特徴とする請求項６に記載の計算装置。
【請求項９】
前記自然数niとして複数の自然数p1，p2，…，pk (k≧2)を用意して、そのそれぞれについて前記写像giを定義し、前記更新回数cを変化させて前記出力部が順次出力するスカラー値の収束速度を取得する収束速度取得部と、
前記収束速度取得部により取得された収束速度が最も速い自然数ph (1≦h≦k)を選択する選択部と、
前記自然数niとして前記自然数phを用いて前記自然数写像giを定義し、前記更新回数cにかえて更新回数c' (c'>c)を用いて前記出力部にスカラー値を出力させる出力制御部と、
をさらに備えることを特徴とする請求項６に記載の計算装置。
【請求項１０】
前記出力部は、前記更新部による更新が行われる度に、前記スカラー値sを計算し、当該スカラー値sと、前回の更新の際におけるこのスカラー値と、を比較して、これが所定の終了条件を満たす場合に、当該スカラー値sを出力する
ことを特徴とする請求項１から９のいずれか１項に記載の計算装置。
【請求項１１】
ベクトル写像fを用いてm (m≧1)次元の非有界領域M = (-∞,+∞)×…×(-∞,+∞)に対してm次元被積分関数Aの多重積分を計算する計算方法であって、
当該計算方法は、第１の記憶部と、第２の記憶部と、第１の計算部と、第２の計算部と、更新部と、出力部と、を備える計算装置にて実行され、
当該ベクトル写像fは、実数の成分からなるm (m≧1)次元のベクトルを実数の成分からなるm次元のベクトルに変換し、
当該ベクトル写像fを任意のｍ次元のベクトルx∈Mに適用した結果f(x)は当該ベクトルxから計算可能であり、
当該ベクトル写像fを所定のm次元のベクトルuに繰り返し適用した結果u，f(u)，f(f(u))，f(f(f(u)))，… の極限分布の多次元密度関数ρを任意のｍ次元のベクトルx∈Mに適用した結果ρ(x)は当該ベクトルxから計算可能であり、
前記第１の記憶部は、m次元のベクトルx=(x1,x2,…,xm)を記憶し、
前記第２の記憶部は、スカラー値wを記憶し、
当該計算方法は、
前記第１の計算部が、前記第１の記憶部に記憶されたベクトルxに、前記ベクトル写像fを適用した結果のベクトルx'=f(x)=(x'1,x'2,…,x'm)を計算する第１の計算工程と、
前記第２の計算部が、前記第１の記憶部に記憶されたベクトルxから、スカラー値w'=A(x)/ρ(x)を計算する第２の計算工程と、
前記更新部が、前記第１の計算工程にて計算されたベクトルx'を前記第１の記憶部に記憶されるベクトルxに代入して、前記第１の記憶部に記憶されるベクトルxを更新し、
前記第２の計算工程にて計算されたスカラー値w'を前記第２の記憶部に記憶されるスカラー値wに加算して、前記第２の記憶部に記憶されるスカラー値wを更新する更新工程と、
前記出力部が、前記更新工程による更新回数がc (c≧1)回となった場合、前記第２の記憶部に記憶されたスカラー値wから、スカラー値s=w/(c+1)を計算して、当該スカラー値sを多重積分の結果として出力する出力工程と、
を備え、
前記第２の記憶部に最初に記憶されるスカラー値は、前記第１の記憶部に最初に記憶されるm次元のベクトルに関数Aを適用したものを、前記第１の記憶部に最初に記憶されるm次元のベクトルに密度関数ρを適用したもので除算した結果である
ことを特徴とする計算方法。
【請求項１２】
ベクトル写像fを用いてm (m≧1)次元の非有界領域M = (-∞,+∞)×…×(-∞,+∞)に対してm次元被積分関数Aの多重積分を計算する計算方法であって、
当該計算方法は、第１の記憶部と、第２の記憶部と、第１の計算部と、第２の計算部と、更新部と、出力部と、を備える計算装置にて実行され、
当該ベクトル写像fは、実数の成分からなるm (m≧1)次元のベクトルを実数の成分からなるm次元のベクトルに変換し、
当該ベクトル写像fを任意のｍ次元のベクトルx∈Mに適用した結果f(x)は当該ベクトルxから計算可能であり、
当該ベクトル写像fを所定のm次元のベクトルuに繰り返し適用した結果u，f(u)，f(f(u))，f(f(f(u)))，… の極限分布の多次元密度関数ρを任意のｍ次元のベクトルx∈Mに適用した結果ρ(x)は当該ベクトルxから計算可能であり、
前記第１の記憶部は、m次元のベクトルx=(x1,x2,…,xm)を記憶し、
前記第２の記憶部は、スカラー値wを記憶し、
当該計算方法は、
前記第１の計算部が、前記第１の記憶部に記憶されたベクトルxに、前記ベクトル写像fを適用した結果のベクトルx'=f(x)=(x'1,x'2,…,x'm)を計算する第１の計算工程と、
前記第２の計算部が、前記第１の記憶部に記憶されたベクトルxから、スカラー値w'=A(x)/ρ(x)を計算する第２の計算工程と、
前記更新部が、前記第１の計算工程にて計算されたベクトルx'を前記第１の記憶部に記憶されるベクトルxに代入して、前記第１の記憶部に記憶されるベクトルxを更新し、
前記第２の計算工程にて計算されたスカラー値w'を前記第２の記憶部に記憶されるスカラー値wに加算して、前記第２の記憶部に記憶されるスカラー値wを更新する更新工程と、
前記出力部が、前記更新工程による更新回数がc (c≧1)回となった場合、前記第２の記憶部に記憶されたスカラー値wから、スカラー値s'=w/c)を計算して、当該スカラー値s'を多重積分の結果として出力する出力工程と、
を備え、
前記第２の記憶部に最初に記憶されるスカラー値は０である
ことを特徴とする計算方法。
【請求項１３】
前記ベクトル写像fは、m次元ベクトルu = (u1，u2，…，um)とm個の写像f1，f2，…，fmとにより、
f(u) = (f1(u1)，f2(u2)，…，fm(um)))
と定義され、当該m個の写像f1，f2，…，fmはいずれもカオス写像である
ことを特徴とする請求項１１または１２に記載の計算方法。
【請求項１４】
前記計算装置は、収束速度取得部と、選択部と、をさらに備え、
当該計算方法は、
前記収束速度取得部が、前記ベクトル写像fとして複数のベクトル写像g1，g2，…，gk (k≧2)を用意し、そのそれぞれについて、前記更新回数cを変化させて前記出力工程にて順次出力されるスカラー値の収束速度を取得する収束速度取得工程と、
前記選択部が、前記収束速度取得工程にて取得された収束速度が最も速いベクトル写像gh (1≦h≦k)を選択する選択工程をさらに備え、
前記ベクトル写像fとして前記ベクトル写像ghを用いて、
前記出力工程は、前記更新回数cにかえて更新回数c' (c'>c)を用いてスカラー値を出力する
ことを特徴とする請求項１１から１３のいずれか１項に記載の計算方法。
【請求項１５】
前記ベクトル写像fを所定のm次元のベクトルu=(u1,u2,…,um)に０回以上適用して得られるベクトル列
u，f(u)，f(f(u))，f(f(f(u)))，…
の極限分布の多次元密度関数ρは、性質
ρ(u) = Πi=0m ρi(ui)；
ρi(ui) 〜 c-i|ui|-(1+a)   for ui → -∞；
ρi(ui) 〜 c+i|ui|-(1+a)   for ui → +∞；
(a>0，1≦i≦m，c-i>0，c+i>0)
を満たす
ことを特徴とする請求項１１から１３のいずれか１項に記載の計算方法。
【請求項１６】
前記ベクトル写像ｆは、1≦i≦mについて定義される関数
fi(t) = gi(di t)/di   (di>0)
により
f(u)=(f1(u1)，f2(u2)，…，fm(um))
と定義され、写像giは、以下の写像φj (1≦j≦8)のいずれかと自然数ni (ni≧2)とにより、
gi(φj(θ)) = φj(niθ)；
φ1(θ) = -sgn(tanθ)/|tanθ|1/a；
φ2(θ) = -sgn(tanθ)×|tanθ|1/a；
φ3(θ) = -sgn(cosθ)/|tanθ|1/a；
φ4(θ) = -sgn(cosθ)×|tanθ|1/a；
φ5(θ) =  sgn(cosθ)/|tanθ|1/a；
φ6(θ) =  sgn(cosθ)×|tanθ|1/a；
φ7(θ) =  sgn(sinθ)/|tanθ|1/a；
φ8(θ) =  sgn(sinθ)×|tanθ|1/a；
sgn(t) =  1   for t>0；
sgn(t) =  0   for t=0；
sgn(t) = -1   for t<0
と定義される
ことを特徴とする請求項１５に記載の計算方法。
【請求項１７】
前記計算装置は、収束速度取得部と、選択部と、出力制御部と、をさらに備え、
当該計算方法は、
前記収束速度取得部が、前記定数aとして複数の正数q1，q2，…，qk (k≧2)を用意して、そのそれぞれについて前記写像fを定義し、前記更新回数cを変化させて前記出力工程にて順次出力されるスカラー値の収束速度を取得する収束速度取得工程と、
前記選択部が、前記収束速度取得工程にて取得された収束速度が最も速い正数qh (1≦h≦k)を選択する選択工程と、
前記出力制御部が、前記定数aとして前記正数qhを用いて前記写像fを定義し、前記更新回数cにかえて更新回数c' (c'>c)を用い、前記出力工程によりスカラー値を出力させる出力制御工程と、
をさらに備えることを特徴とする請求項１５または１６に記載の計算方法。
【請求項１８】
前記計算装置は、収束速度取得部と、選択部と、出力制御部と、をさらに備え、
当該計算方法は、
前記収束速度取得部が、前記写像giを前記写像φjのいずれか複数を用いて定義し、そのそれぞれについて、前記更新回数cを変化させて前記出力工程にて順次出力されるスカラー値の収束速度を取得する収束速度取得工程と、
前記選択部が、前記収束速度取得工程にて取得された収束速度が最も速い写像φjを選択する選択工程と、
前記出力制御部が、前記写像giを前記選択工程にて選択された写像φjを用いて定義し、前記更新回数cにかえて更新回数c' (c'>c)を用いて前記出力工程によりスカラー値を出力させる出力制御工程と、
をさらに備えることを特徴とする請求項１６に記載の計算方法。
【請求項１９】
前記計算装置は、収束速度取得部と、選択部と、出力制御部と、をさらに備え、
当該計算方法は、
前記収束速度取得部が、前記自然数niとして複数の自然数p1，p2，…，pk (k≧2)を用意して、そのそれぞれについて前記写像giを定義し、前記更新回数cを変化させて前記出力工程にて順次出力されるスカラー値の収束速度を取得する収束速度取得工程と、
前記選択部が、前記収束速度取得工程にて取得された収束速度が最も速い自然数ph (1≦h≦k)を選択する選択工程と、
前記出力制御部が、前記自然数niとして前記自然数phを用いて前記自然数写像giを定義し、前記更新回数cにかえて更新回数c' (c'>c)を用いて前記出力工程によりスカラー値を出力させる出力制御工程と、
をさらに備えることを特徴とする請求項１６に記載の計算方法。
【請求項２０】
前記出力工程は、前記更新工程による更新が行われる度に、前記スカラー値sを計算し、当該スカラー値sと、前回の更新の際におけるこのスカラー値と、を比較して、これが所定の終了条件を満たす場合に、当該スカラー値sを出力する
ことを特徴とする請求項１１から１９のいずれか１項に記載の計算方法。
【請求項２１】
ベクトル写像fを用いてm (m≧1)次元の非有界領域M = (-∞,+∞)×…×(-∞,+∞)に対してm次元被積分関数Aの多重積分を計算するためのプログラムを記録した、コンピュータ読取可能な情報記録媒体であって、
当該ベクトル写像fは、実数の成分からなるm (m≧1)次元のベクトルを実数の成分からなるm次元のベクトルに変換し、
当該ベクトル写像fを任意のｍ次元のベクトルx∈Mに適用した結果f(x)は当該ベクトルxから計算可能であり、
当該ベクトル写像fを所定のm次元のベクトルuに繰り返し適用した結果u，f(u)，f(f(u))，f(f(f(u)))，… の極限分布の多次元密度関数ρを任意のｍ次元のベクトルx∈Mに適用した結果ρ(x)は当該ベクトルxから計算可能であり、
当該プログラムは、コンピュータを、
m次元のベクトルx=(x1,x2,…,xm)を記憶する第１の記憶部、
スカラー値wを記憶する第２の記憶部、
前記第１の記憶部に記憶されたベクトルxに、前記ベクトル写像fを適用した結果のベクトルx'=f(x)=(x'1,x'2,…,x'm)を計算する第１の計算部、
前記第１の記憶部に記憶されたベクトルxから、スカラー値w'=A(x)/ρ(x)を計算する第２の計算部、
前記第１の計算部により計算されたベクトルx'を前記第１の記憶部に記憶されるベクトルxに代入して、前記第１の記憶部に記憶されるベクトルxを更新し、
前記第２の計算部により計算されたスカラー値w'を前記第２の記憶部に記憶されるスカラー値wに加算して、前記第２の記憶部に記憶されるスカラー値wを更新する更新部、および、
前記更新部による更新回数がc (c≧1)回となった場合、前記第２の記憶部に記憶されたスカラー値wから、スカラー値s=w/(c+1)を計算して、当該スカラー値sを多重積分の結果として出力する出力部、
として機能させ、
前記第２の記憶部に最初に記憶されるスカラー値は、前記第１の記憶部に最初に記憶されるm次元のベクトルに関数Aを適用したものを、前記第１の記憶部に最初に記憶されるm次元のベクトルに密度関数ρを適用したもので除算した結果である
ように機能させることを特徴とする情報記録媒体。
【請求項２２】
ベクトル写像fを用いてm (m≧1)次元の非有界領域M = (-∞,+∞)×…×(-∞,+∞)に対してm次元被積分関数Aの多重積分を計算するためのプログラムを記録した、コンピュータ読取可能な情報記録媒体であって、
当該ベクトル写像fは、実数の成分からなるm (m≧1)次元のベクトルを実数の成分からなるm次元のベクトルに変換し、
当該ベクトル写像fを任意のｍ次元のベクトルx∈Mに適用した結果f(x)は当該ベクトルxから計算可能であり、
当該ベクトル写像fを所定のm次元のベクトルuに繰り返し適用した結果u，f(u)，f(f(u))，f(f(f(u)))，… の極限分布の多次元密度関数ρを任意のｍ次元のベクトルx∈Mに適用した結果ρ(x)は当該ベクトルxから計算可能であり、
当該プログラムは、コンピュータを、
m次元のベクトルx=(x1,x2,…,xm)を記憶する第１の記憶部、
スカラー値wを記憶する第２の記憶部、
前記第１の記憶部に記憶されたベクトルxに、前記ベクトル写像fを適用した結果のベクトルx'=f(x)=(x'1,x'2,…,x'm)を計算する第１の計算部、
前記第１の記憶部に記憶されたベクトルxから、スカラー値w'=A(x)/ρ(x)を計算する第２の計算部、
前記第１の計算部により計算されたベクトルx'を前記第１の記憶部に記憶されるベクトルxに代入して、前記第１の記憶部に記憶されるベクトルxを更新し、
前記第２の計算部により計算されたスカラー値w'を前記第２の記憶部に記憶されるスカラー値wに加算して、前記第２の記憶部に記憶されるスカラー値wを更新する更新部、および、
前記更新部による更新回数がc (c≧1)回となった場合、前記第２の記憶部に記憶されたスカラー値wから、スカラー値s'=w/cを計算して、当該スカラー値s'を多重積分の結果として出力する出力部、
として機能させ、
前記第２の記憶部に最初に記憶されるスカラー値は０であるように機能させることを特徴とする情報記録媒体。
【請求項２３】
前記プログラムは、前記コンピュータにおいて、
前記ベクトル写像fは、m次元ベクトルu = (u1，u2，…，um)とm個の写像f1，f2，…，fmとにより、
f(u) = (f1(u1)，f2(u2)，…，fm(um)))
と定義され、当該m個の写像f1，f2，…，fmはいずれもカオス写像である
ように機能させることを特徴とする請求項２１または２２に記載の情報記録媒体。
【請求項２４】
前記プログラムは、前記コンピュータを、
前記ベクトル写像fとして複数のベクトル写像g1，g2，…，gk (k≧2)を用意し、そのそれぞれについて、前記更新回数cを変化させて前記出力部が順次出力するスカラー値の収束速度を取得する収束速度取得部、
前記収束速度取得部により取得された収束速度が最も速いベクトル写像gh (1≦h≦k)を選択する選択部、および、
前記ベクトル写像fとして前記ベクトル写像ghを用いて、前記更新回数cにかえて更新回数c' (c'>c)を用いて前記出力部にスカラー値を出力させる出力制御部、
としてさらに機能させることを特徴とする請求項２１から２３のいずれか１項に記載の情報記録媒体。
【請求項２５】
前記ベクトル写像fを所定のm次元のベクトルu=(u1,u2,…,um)に０回以上適用して得られるベクトル列
u，f(u)，f(f(u))，f(f(f(u)))，…
の極限分布の多次元密度関数ρは、性質
ρ(u) = Πi=0m ρi(ui)；
ρi(ui) 〜 c-i|ui|-(1+a)   for ui → -∞；
ρi(ui) 〜 c+i|ui|-(1+a)   for ui → +∞；
(a>0，1≦i≦m，c-i>0，c+i>0)
を満たす
ことを特徴とする請求項２１から２３のいずれか１項に記載の情報記録媒体。
【請求項２６】
前記ベクトル写像ｆは、1≦i≦mについて定義される関数
fi(t) = gi(di t)/di   (di>0)
により
f(u)=(f1(u1)，f2(u2)，…，fm(um))
と定義され、写像giは、以下の写像φj (1≦j≦8)のいずれかと自然数ni (ni≧2)とにより、
gi(φj(θ)) = φj(niθ)；
φ1(θ) = -sgn(tanθ)/|tanθ|1/a；
φ2(θ) = -sgn(tanθ)×|tanθ|1/a；
φ3(θ) = -sgn(cosθ)/|tanθ|1/a；
φ4(θ) = -sgn(cosθ)×|tanθ|1/a；
φ5(θ) =  sgn(cosθ)/|tanθ|1/a；
φ6(θ) =  sgn(cosθ)×|tanθ|1/a；
φ7(θ) =  sgn(sinθ)/|tanθ|1/a；
φ8(θ) =  sgn(sinθ)×|tanθ|1/a；
sgn(t) =  1   for t>0；
sgn(t) =  0   for t=0；
sgn(t) = -1   for t<0
と定義される
ことを特徴とする請求項２５に記載の情報記録媒体。
【請求項２７】
前記プログラムは、前記コンピュータを、
前記定数aとして複数の正数q1，q2，…，qk (k≧2)を用意して、そのそれぞれについて前記写像fを定義し、前記更新回数cを変化させて前記出力部が順次出力するスカラー値の収束速度を取得する収束速度取得部、
前記収束速度取得部により取得された収束速度が最も速い正数qh (1≦h≦k)を選択する選択部、および、
前記定数aとして前記正数qhを用いて前記写像fを定義し、前記更新回数cにかえて更新回数c' (c'>c)を用いて前記出力部にスカラー値を出力させる出力制御部、
としてさらに機能させることを特徴とする請求項２５または２６に記載の情報記録媒体。
【請求項２８】
前記プログラムは、前記コンピュータを、
前記写像giを前記写像φjのいずれか複数を用いて定義し、そのそれぞれについて、前記更新回数cを変化させて前記出力部が順次出力するスカラー値の収束速度を取得する収束速度取得部、
前記収束速度取得部により取得された収束速度が最も速い写像φjを選択する選択部、および、
前記写像giを前記選択部により選択された写像φjを用いて定義し、前記更新回数cにかえて更新回数c' (c'>c)を用いて前記出力部にスカラー値を出力させる出力制御部、
としてさらに機能させることを特徴とする請求項２６に記載の情報記録媒体。
【請求項２９】
前記プログラムは、前記コンピュータを、
前記自然数niとして複数の自然数p1，p2，…，pk (k≧2)を用意して、そのそれぞれについて前記写像giを定義し、前記更新回数cを変化させて前記出力部が順次出力するスカラー値の収束速度を取得する収束速度取得部、
前記収束速度取得部により取得された収束速度が最も速い自然数ph (1≦h≦k)を選択する選択部、および、
前記自然数niとして前記自然数phを用いて前記自然数写像giを定義し、前記更新回数cにかえて更新回数c' (c'>c)を用いて前記出力部にスカラー値を出力させる出力制御部、
としてさらに機能させることを特徴とする請求項２６に記載の情報記録媒体。
【請求項３０】
前記出力部は、前記更新部による更新が行われる度に、前記スカラー値sを計算し、当該スカラー値sと、前回の更新の際におけるこのスカラー値と、を比較して、これが所定の終了条件を満たす場合に、当該スカラー値sを出力する
ことを特徴とする請求項２１から２９のいずれか１項に記載の情報記録媒体。
【請求項３１】
前記情報記録媒体は、コンパクトディスク、フレキシブルディスク、ハードディスク、光磁気ディスク、ディジタルビデオディスク、磁気テープ、または、半導体メモリである
ことを特徴とする請求項２１から３０のいずれか１項に記載の情報記録媒体。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
本発明は、多重積分の計算装置、計算方法、および、情報記録媒体に関する。
特に、多次元の非有界領域上で定義される関数の積分の近似値を高速に計算する多重積分の計算装置、計算方法、および、これらを実現するプログラムを記録した情報記録媒体に関する。
【０００２】
【従来の技術】
従来から、多次元の非有界領域上の多重積分を求める問題は、金融デリバティブ計算のリスク評価計算や、ＣＤＭＡ（Code Division Multiple Access；符号分割多重アクセス）などの多重通信の干渉ノイズの評価など、多岐に渡る分野で現れている。
【０００３】
積分の領域が１次元や２次元の場合、これを数値積分するには、積分区間を細かいメッシュに分け、メッシュ内の被積分関数の値をある代表値で置き換え、メッシュの代表値の平均値を計算することにより、近似的な数値積分値を求めることができた。
【０００４】
ところが、積分領域の次元が高くなると、メッシュの数が次元数mに関しm乗のオーダーで指数関数的に増える問題があった。この問題を「次元の呪（curse of dimensionality）」と呼ぶ。
【０００５】
この問題を解決する手法の一つとして、モンテカルロ計算による数値積分法が提案されている。モンテカルロ計算では、擬似乱数（low-discrepancy sequences等の準乱数を含む。）を用いるが、この擬似乱数は、できるだけ一様に近い分布をすることが望ましい。
【０００６】
しかしながら、積分領域が無限区間(-∞,+∞)などの非有界領域の場合には、有限区間上の一様乱数を当該非有界領域上へ線形変換して正規化することはできない。そこで、何らかの非線形変換をする必要がある。
【０００７】
J. von NeumannとS. Ulamは、１９４７年に、一様乱数から任意の所望の密度関数に従う乱数を発生させる手法を発明し、核分裂における各中性子の軌跡のコンピュータによる模擬実験に成功した。この手法は、一般に「逆関数法」と呼ばれる。密度関数ρの積分関数を
η(X) = ∫C ρ(x)dx   (C = {u | -∞≦u≦X})
とするとき、一様乱数への変換に逆関数η-1(X)を用いるからである。
【０００８】
一方で、近年のカオス理論の発展により、有限区間上の非一様な分布を持つ乱数を得ることができることが判明しており、乱数発生における写像力学系の重要性はますます高まってきている。
【０００９】
このような乱数発生には、三角関数を含む楕円関数の加法定理から導かれる有理写像を漸化式に用いており、以下のような有利な性質があることがわかっている。
【００１０】
すなわち、出力される乱数列には、わずかな例外を除き、周期がないため、同じ列が繰り返し出力されることがない。
【００１１】
また、カオスの持つ初期値鋭敏性により、乱数の種（漸化式に与える初期値）をわずかに変えるだけで、得られる乱数列は全く違ったものになる。
【００１２】
さらに、乱数の分布を表す密度関数が既知の解析的関数である。
【００１３】
このような写像としては、ウラム＝フォン・ノイマン写像
f(x) = 4x(1-x)
キュービック写像
f(x) = x(3-4x)2
クインティック写像
f(x) = x(5-20x+16x2)2
などが知られている。
【００１４】
これらの写像のいずれを選んだ場合であっても、適当な初期値ξ (0＜ξ＜1)を与え、漸化式
x[0] = ξ；
x[i+1] = f(x[i])  (i≧0)
により乱数列x[0], x[1], x[2], …を得た場合、この乱数列x[i] (i≧0)の極限分布を表す密度関数は、
ρ(x)=1/(π(x(1-x))1/2)
となる。
【００１５】
これらの写像を用いた有界領域上の乱数を生成する手法と、そのモンテカルロ数値積分への応用については、本願の発明者らによる出願に係る特開平１０−２８３３４４号公報に開示されている。また、その理論的背景については以下の文献に開示されている。
S. M. Ulam and J. von Neumann, Bull. Math. Soc. 53 (1947) p1120.
R. L. Adler and T. J. Rivlin, Proc. Am. Math. Soc. 15 (1964) p794.
K. Umeno, Method of constructing exactly solvable chaos, Phys. Rev. E (1997) Vol.55, pp5280-5284.
【００１６】
【発明が解決しようとする課題】
しかしながら、上記のような従来の手法には、以下のような問題があった。
【００１７】
すなわち、非有界領域上におけるモンテカルロ数値積分法では、非有界領域上に既知の解析関数で与えられる密度関数に従って分布する乱数が必要である。従って一様乱数を得た後、必ず、非有界領域上に分布する乱数に変換をしなければならない。上述のJ. von Neumannの逆関数法により、所望の密度関数に分布する乱数を得るには、その密度関数の積分関数を求めてから、更にその積分関数の逆関数を求めなければならない。
【００１８】
しかしながら、そもそも密度関数の積分を解析的に求めるのが、ほとんどの場合不可能である。また、その密度関数の積分関数の逆関数を求める際に、数値計算に伴う近似と計算時間を考慮しなければならない。モンテカルロ計算では、乱数を大量に発生させる必要があるため、この手法で積分を求めることは難しい。
【００１９】
また、通常の株価からオプションの値を出すブラック・ショールズ方程式を解く場合にも、非有界領域における数値積分が必要である。株価変動が正規分布（ガウス分布）に従うと仮定すれば、非有解領域上のガウス積分を用いて解を得ることができる。
【００２０】
しかしながら、株価変動を実際に観測した結果、その分布の密度関数は、羃的な裾野を持つことがわかっている（R. N. Mantegna and H. E. Stanley, Nature, vol. 376, pp.46-49, 1995）。羃的裾野を持つ密度関数の積分は、ガウス分布とは異なり、初等関数の形で書けない。このため、既存の一様乱数を基礎とする非有界モンテカルロ計算法では、解を得ることが極めて難しい。
【００２１】
この他、冪的な待ち時間分布を持つ通信トラフィック、ＣＤＭＡにおける干渉雑音の評価にも、同様の問題が生じていた。
【００２２】
したがって、非有界領域上の多重積分を効率良く高速に計算する計算手法に対する要望は極めて大きく、その実現は産業界から広く望まれている。
【００２３】
本発明は、多次元の非有界領域上で定義される関数の積分の近似値を高速に計算するのに好適な、多重積分の計算装置、計算方法、および、これらを実現するプログラムを記録した情報記録媒体を提供することを目的とする。
【００２４】
【課題を解決するための手段】
以上の目的を達成するため、本発明の原理にしたがって、下記の発明を開示する。
【００２５】
本発明の第１の観点に係る多重積分の計算装置は、ベクトル写像fを用いて、多次元被積分関数Aの多重積分を計算する。ここで、ベクトル写像fは、実数の成分からなるm (m≧1)次元のベクトルを実数の成分からなるm次元のベクトルに変換し、非有界領域のサポートを持ち、所定のm次元のベクトルuに写像fを繰り返し適用した結果の極限分布の多次元密度関数ρは解析的に得られる。
【００２６】
図１は、本計算装置の概要構成を示す模式図である。以下、本図を参照して説明する。
【００２７】
本計算装置１０１は、第１の記憶部１１１と、第２の記憶部１１２と、第１の計算部１２１と、第２の計算部１２２と、更新部１３１と、出力部１４１と、を備える。
【００２８】
ここで、第１の記憶部１１１は、m次元のベクトルx=(x1,x2,…,xm)を記憶する。
【００２９】
一方、第２の記憶部１１２は、スカラー値wを記憶する。
【００３０】
そして、第１の計算部１２１は、第１の記憶部１１１に記憶されたベクトルxに、ベクトル写像fを適用した結果のベクトルx'=f(x)=(x'1,x'2,…,x'm)を計算する。
【００３１】
さらに、第２の計算部１２２は、第１の記憶部１１１に記憶されたベクトルxと、第２の記憶部１１２に記憶されたスカラー値wとから、スカラー値w'=A(x)/ρ(x)を計算する。
【００３２】
一方、更新部１３１は、第１の計算部１２１により計算されたベクトルx'を第１の記憶部１１１に記憶させて更新し、第２の計算部１２２により計算されたスカラー値w'を第２の記憶部１１２に記憶される値に加算して、これを更新する。
【００３３】
そして、出力部１４１は、更新部１３１による更新回数がc (c≧1)回となった場合、第２の記憶部１１２に記憶されたスカラー値wから、スカラー値s=w/(c+1)を計算して、当該スカラー値sを多重積分の結果として出力する。
【００３４】
本発明は、多次元の非有界領域Mについて、エルゴード性
∫M A(x)dx
= ∫M (A(x)/ρ(x))・ρ(x)dx
= (1/N)・Σj=1N A(X[j])/ρ(X[j])  (N → ∞)；
X[j+1] = f(X[j])  (j≧0)；
X[0]は適当なm次元ベクトル
が成立することによる。
【００３５】
この場合、計算誤差は、最悪時であっても、次元数mに依存せず(1/N)1/2のオーダーである。また、発明者による最近の研究成果（K. Umeno「Chaotic Monte Calro Computation： a Dynamical Effect of Random-Number Generations」to appear in Japanese Journal of Applied Physics, March 2000）より、Super-efficiencyと呼ばれる条件を満たす場合は、次元数mに依存せず1/Nのオーダーとなることが判明している。したがって、高速な積分計算が可能である。
【００３６】
また、本発明の計算装置は、第２の記憶部に最初に記憶されるスカラー値は、第１の記憶部に最初に記憶されるｍ次元のベクトルに関数Aを適用したものを、第１の記憶部に最初に記憶されるｍ次元のベクトルに密度関数ρを適用したもので除算した結果であるように構成することができる。
【００３７】
また、本発明の計算装置において、第２の記憶部に最初に記憶されるスカラー値は０であり、出力部は、スカラー値sにかえて、スカラー値s'=w/cを計算して、当該スカラー値s'を多重積分の結果として出力するように構成することができる。
【００３８】
また、本発明の計算装置は、収束速度取得部と、選択部と、出力制御部と、をさらに備えるように構成することができる。
【００３９】
ここで、収束速度取得部は、ベクトル写像fとして複数のベクトル写像g1, g2,…, gk (k≧2)を用意し、そのそれぞれについて、更新回数cを変化させて出力部が順次出力するスカラー値の収束速度を取得する。
【００４０】
一方、選択部は、収束速度取得部により取得された収束速度が最も速いベクトル写像gh (1≦h≦k)を選択する。
【００４１】
さらに、出力制御部は、ベクトル写像fとしてベクトル写像ghを用いて、更新回数cにかえて更新回数c' (c'＞c)を用いて出力部にスカラー値を出力させる。
【００４２】
また、本発明の計算装置において、ベクトル写像fを所定のm次元のベクトルu=(u1,u2,…,um)に０回以上適用して得られるベクトル列
u, f(u), f(f(u)), f(f(f(u))), …
の極限分布の多次元密度関数ρは、性質
ρ(u) = Πi=0m ρi(ui)；
ρi(ui) 〜 c-i|ui|-(1+a)   for ui → -∞；
ρi(ui) 〜 c+i|ui|-(1+a)   for ui → +∞；
(a＞0, 1≦i≦m, c-i＞0, c+i＞0)
を満たすように構成することができる。
【００４３】
また、本発明の計算装置において、ベクトル写像ｆは、1≦i≦mについて定義される関数
fi(t) = gi(di t)/di   (di＞0)
により
f(u)=(f1(u1), f2(u2), …, fm(um))
と定義され、写像giは、以下の写像φj (1≦j≦8)のいずれかと自然数ni (ni≧2)とにより、
gi(φj(θ)) = φj(niθ)；
φ1(θ) = -sgn(tanθ)/|tanθ|1/a；
φ2(θ) = -sgn(tanθ)×|tanθ|1/a；
φ3(θ) = -sgn(cosθ)/|tanθ|1/a；
φ4(θ) = -sgn(cosθ)×|tanθ|1/a；
φ5(θ) =  sgn(cosθ)/|tanθ|1/a；
φ6(θ) =  sgn(cosθ)×|tanθ|1/a；
φ7(θ) =  sgn(sinθ)/|tanθ|1/a；
φ8(θ) =  sgn(sinθ)×|tanθ|1/a；
sgn(t) =  1   for t＞0；
sgn(t) =  0   for t=0；
sgn(t) = -1   for t＜0
と定義されるように構成することができる。
【００４４】
また、本発明の計算装置は、収束速度取得部と、選択部と、出力制御部と、を備えるように構成することができる。
【００４５】
ここで、収束速度取得部は、定数aとして複数の正数q1, q2, …, qk (k≧2)を用意して、そのそれぞれについて写像fを定義し、更新回数cを変化させて出力部が順次出力するスカラー値の収束速度を取得する。
【００４６】
一方、選択部は、収束速度取得部により取得された収束速度が最も速い正数qh (1≦h≦k)を選択する。
【００４７】
さらに、出力制御部は、定数aとして正数qhを用いて写像fを定義し、更新回数cにかえて更新回数c' (c'＞c)を用いて出力部にスカラー値を出力させる。
【００４８】
また、本発明の計算装置は、収束速度取得部と、選択部と、出力制御部と、を備えるように構成することができる。
【００４９】
ここで、収束速度取得部は、写像giを写像φjのいずれか複数を用いて定義し、そのそれぞれについて、更新回数cを変化させて出力部が順次出力するスカラー値の収束速度を取得する。
【００５０】
一方、選択部は、収束速度取得部により取得された収束速度が最も速い写像φjを選択する。
【００５１】
さらに、出力制御部は、写像giを選択部により選択された写像φjを用いて定義し、更新回数cにかえて更新回数c' (c'＞c)を用いて出力部にスカラー値を出力させる。
【００５２】
また、本発明の計算装置は、収束速度取得部と、選択部と、出力制御部と、を備えるように構成することができる。
【００５３】
ここで、収束速度取得部は、自然数niとして複数の自然数p1, p2, …, pk (k≧2)を用意して、そのそれぞれについて写像giを定義し、更新回数cを変化させて出力部が順次出力するスカラー値の収束速度を取得する。
【００５４】
一方、選択部は、収束速度取得部により取得された収束速度が最も速い自然数ph (1≦h≦k)を選択する。
【００５５】
さらに、出力制御部は、自然数niとして自然数phを用いて自然数写像giを定義し、更新回数cにかえて更新回数c' (c'＞c)を用いて出力部にスカラー値を出力させる。
【００５６】
また、本発明の計算装置の出力部は、更新部による更新が行われる度に、スカラー値sを計算し、当該スカラー値sと、前回の更新の際におけるこのスカラー値と、を比較して、これが所定の終了条件を満たす場合に、当該スカラー値sを出力するように構成することができる。
【００５７】
本発明の第１の観点に係る多重積分の計算方法は、ベクトル写像fと、m次元のベクトルx=(x1,x2,…,xm)を記憶する第１の記憶部と、スカラー値wを記憶する第２の記憶部と、を用いて、多次元被積分関数Aの多重積分を計算する。
【００５８】
ここで、ベクトル写像fは、実数の成分からなるm (m≧1)次元のベクトルを実数の成分からなるm次元のベクトルに変換し、非有界領域のサポートを持ち、所定のm次元のベクトルuに写像fを繰り返し適用した結果の極限分布の多次元密度関数ρは解析的に得られる。
【００５９】
本計算方法は、第１の計算工程と、第２の計算工程と、更新工程と、出力工程と、を備えるように構成する。
【００６０】
ここで、第１の計算工程では、第１の記憶部に記憶されたベクトルxに、ベクトル写像fを適用した結果のベクトルx'=f(x)=(x'1,x'2,…,x'm)を計算する。
【００６１】
一方、第２の計算工程では、第１の記憶部に記憶されたベクトルxと、第２の計算部に記憶されたスカラー値wとから、スカラー値w'=A(x)/ρ(x)を計算する。
【００６２】
さらに、更新工程では、第１の計算工程にて計算されたベクトルx'を第１の記憶部に記憶させて更新し、第２の計算工程にて計算されたスカラー値w'を第２の記憶部に記憶される値に加算してこれを更新する。
【００６３】
そして、出力工程では、更新工程による更新回数がc (c≧1)回となった場合、第２の記憶部に記憶されたスカラー値wから、スカラー値s=w/(c+1)を計算して、当該スカラー値sを多重積分の結果として出力する。
【００６４】
また、本発明の計算方法において、第２の記憶部に最初に記憶されるスカラー値は、第１の記憶部に最初に記憶されるｍ次元のベクトルに関数Aを適用したものを、第１の記憶部に最初に記憶されるｍ次元のベクトルに密度関数ρを適用したもので除算した結果であるように構成することができる。
【００６５】
また、本発明の計算方法において、第２の記憶部に最初に記憶されるスカラー値は０であり、出力工程は、スカラー値sにかえて、スカラー値s'=w/cを計算して、当該スカラー値s'を多重積分の結果として出力するように構成することができる。
【００６６】
また、本発明の計算方法は、収束速度取得工程と、選択工程とをさらに備えるように構成することができる。
【００６７】
ここで、収束速度取得工程では、ベクトル写像fとして複数のベクトル写像g1, g2, …, gk (k≧2)を用意し、そのそれぞれについて、更新回数cを変化させて出力工程にて順次出力されるスカラー値の収束速度を取得する。
【００６８】
一方、選択工程では、収束速度取得工程にて取得された収束速度が最も速いベクトル写像gh (1≦h≦k)を選択する。
【００６９】
また、ベクトル写像fとしてベクトル写像ghを用いる。
【００７０】
さらに、出力工程では、更新回数cにかえて更新回数c' (c'＞c)を用いてスカラー値を出力する。
【００７１】
また、本発明の計算方法において、ベクトル写像fを所定のm次元のベクトルu=(u1,u2,…,um)に０回以上適用して得られるベクトル列
u, f(u), f(f(u)), f(f(f(u))), …
の極限分布の多次元密度関数ρは、性質
ρ(u) = Πi=0m ρi(ui)；
ρi(ui) 〜 c-i|ui|-(1+a)   for ui → -∞；
ρi(ui) 〜 c+i|ui|-(1+a)   for ui → +∞；
(a＞0, 1≦i≦m, c-i＞0, c+i＞0)
を満たすように構成することができる。
【００７２】
また、本発明の計算方法において、ベクトル写像ｆは、1≦i≦mについて定義される関数
fi(t) = gi(di t)/di   (di＞0)
により
f(u)=(f1(u1), f2(u2), …, fm(um))
と定義され、写像giは、以下の写像φj (1≦j≦8)のいずれかと自然数ni (ni≧2)とにより、
gi(φj(θ)) = φj(niθ)；
φ1(θ) = -sgn(tanθ)/|tanθ|1/a；
φ2(θ) = -sgn(tanθ)×|tanθ|1/a；
φ3(θ) = -sgn(cosθ)/|tanθ|1/a；
φ4(θ) = -sgn(cosθ)×|tanθ|1/a；
φ5(θ) =  sgn(cosθ)/|tanθ|1/a；
φ6(θ) =  sgn(cosθ)×|tanθ|1/a；
φ7(θ) =  sgn(sinθ)/|tanθ|1/a；
φ8(θ) =  sgn(sinθ)×|tanθ|1/a；
sgn(t) =  1   for t＞0；
sgn(t) =  0   for t=0；
sgn(t) = -1   for t＜0
と定義されるように構成することができる。
【００７３】
また、本発明の計算方法は、収束速度取得工程と、選択工程と、出力制御工程と、をさらに備えるように構成することができる。
【００７４】
ここで、収束速度取得工程では、定数aとして複数の正数q1, q2, …, qk (k≧2)を用意して、そのそれぞれについて写像fを定義し、更新回数cを変化させて出力工程にて順次出力されるスカラー値の収束速度を取得する。
【００７５】
一方、選択工程では、収束速度取得工程にて取得された収束速度が最も速い正数qh (1≦h≦k)を選択する。
【００７６】
さらに、出力制御工程では、定数aとして正数qhを用いて写像fを定義し、更新回数cにかえて更新回数c' (c'＞c)を用い、出力工程によりスカラー値を出力させる。
【００７７】
また、本発明の計算方法は、収束速度取得工程と、選択工程と、出力制御工程と、をさらに備えるように構成することができる。
【００７８】
ここで、収束速度取得工程では、写像giを写像φjのいずれか複数を用いて定義し、そのそれぞれについて、更新回数cを変化させて出力工程にて順次出力されるスカラー値の収束速度を取得する。
【００７９】
一方、選択工程では、収束速度取得工程にて取得された収束速度が最も速い写像φjを選択する。
【００８０】
さらに、出力制御工程では、写像giを選択工程にて選択された写像φjを用いて定義し、更新回数cにかえて更新回数c' (c'＞c)を用いて出力工程によりスカラー値を出力させる。
【００８１】
また、本発明の計算方法は、収束速度取得工程と、選択工程と、出力制御工程と、をさらに備えるように構成することができる。
【００８２】
ここで、収束速度取得工程では、自然数niとして複数の自然数p1, p2, …, pk (k≧2)を用意して、そのそれぞれについて写像giを定義し、更新回数cを変化させて出力工程にて順次出力されるスカラー値の収束速度を取得する。
【００８３】
一方、選択工程では、収束速度取得工程にて取得された収束速度が最も速い自然数ph (1≦h≦k)を選択する。
【００８４】
さらに、出力制御工程では、自然数niとして自然数phを用いて自然数写像giを定義し、更新回数cにかえて更新回数c' (c'＞c)を用いて出力工程によりスカラー値を出力させる。
【００８５】
また、本発明の計算方法において、出力工程は、更新工程による更新が行われる度に、スカラー値sを計算し、当該スカラー値sと、前回の更新の際におけるこのスカラー値と、を比較して、これが所定の終了条件を満たす場合に、当該スカラー値sを出力するように構成することができる。
【００８６】
本発明の多重積分の計算装置と、計算方法とを実現するプログラムをコンパクトディスク、フレキシブルディスク、ハードディスク、光磁気ディスク、ディジタルビデオディスク、磁気テープ、半導体メモリなどの情報記録媒体に記録することができる。
【００８７】
本発明の情報記録媒体に記録されたプログラムを、記憶装置、計算装置、出力装置などを備える情報処理装置、たとえば汎用コンピュータ、ゲーム装置、携帯情報端末、移動体電話で実行することにより、上記の多重積分の計算装置と、計算方法とを実現することができる。
【００８８】
また、情報処理装置とは独立して、本発明のプログラムを記録した情報記録媒体を配布、販売することができる。
【００８９】
【発明の実施の形態】
以下に本発明の実施例を説明する。なお、以下に説明する実施例は説明のためのものであり、本願発明の範囲を制限するものではない。したがって、当業者であればこれらの各要素もしくは全要素をこれと均等なものに置換した実施形態を採用することが可能であり、これらの実施形態も本願発明の範囲に含まれる。
【００９０】
（第１の実施形態）
図２は、本発明の多重積分の計算装置を、汎用コンピュータなどの情報処理装置において実現する際の、当該情報処理装置の概要構成を示す模式図である。以下、本図を参照して説明する。
【００９１】
情報処理装置３０１は、ＣＰＵ（Central Processing Unit；中央処理ユニット）３０２により制御される。ＲＡＭ（Random Access Memory；ランダム・アクセス・メモリ）などの主記憶装置３０３には、一時的なデータなどが記憶される。ハードディスク、フロッピーディスク、ＣＤ−ＲＯＭ（Compact Disk Read Only Memory）、磁気テープ、光磁気ディスクなどの外部記憶装置３０４には、ＣＰＵ  ３０２が実行するプログラムが記憶される。
【００９２】
情報処理装置３０１に電源が投入されると、ＣＰＵ  ３０２は、まず、ＲＯＭ（Read Only Memory；読出専用メモリ）３０８に記憶されている初期プログラムローダ（Initial Program Loader；ＩＰＬ）と呼ばれるプログラムを実行する。しかる後に外部記憶装置３０４などからオペレーティングシステムのプログラムやアプリケーションのプログラムなどを主記憶装置３０３にロードして実行する。
【００９３】
実行した結果は、外部記憶装置３０４にファイルとして記憶したり、ＣＲＴ（Cathode Ray Tube）や液晶ディスプレイなどの表示装置３０５に表示することができる。情報処理装置３０１のユーザは、マウスやキーボードなどの入力装置３０６を用いて情報処理装置３０１に対する指示を与える。
【００９４】
ここで、主記憶装置３０３は、第１の記憶部および第２の記憶部として機能する。
【００９５】
一方、ＣＰＵ  ３０２が備えるＡＬＵ（Arithmetic Logic Unit；演算論理ユニット）や、数値演算コプロセッサ（図示せず）は、第１の計算部および第２の計算部として機能する。
【００９６】
さらに、ＣＰＵ  ３０２は、主記憶装置３０３と共働して更新部、および、出力部として機能する。
【００９７】
また、外部記憶装置３０４は、本発明のプログラムを記録する情報記録媒体として機能する。
【００９８】
なお、これらの機能は、下記のように数値演算とメモリの更新により実現されるので、汎用コンピュータなどの情報処理装置を用いるのではなく、専用の電子回路を用いて構成することも容易であり、このような実施形態も本発明の範囲に含まれる。
【００９９】
（計算処理の流れ）
図３は、本発明の計算装置において行われる処理の手順、すなわち、本発明の計算方法の手順の第１の実施例の流れを示すフローチャートである。以下、本図を参照して説明する。なお、本処理は、段落００３４で既に述べたΣを含む式の計算を十分に大きなN（以下では定数Cに対応する。）について行うことに相応するものである。
【０１００】
主記憶装置３０３（典型的にはＲＡＭ）内には、m次元ベクトル変数x領域（第１の記憶部）、スカラー変数w領域（第２の記憶部）、一時的な記憶領域として用いるためのm次元ベクトル変数x'領域とスカラー変数w'領域、更新回数c領域、および、積分の結果変数s領域があらかじめ確保される。以下、これらの領域をそれぞれ、変数x、変数w、変数x'、変数w'、変数c、変数sと略記し、これらにその時に記憶されている値を、それぞれ、x、w、x'、w'、c、sのように略記する。
【０１０１】
また、非有界領域M = (-∞,+∞)×…×(-∞,+∞)に対して積分を行うべき多次元関数をAとし、後述する条件を満たすベクトルカオス写像をf、これを繰り返し適用したその極限分布の密度関数をρとする。
【０１０２】
まず、ＣＰＵ  ３０２は、変数cに値０を記憶させる（ステップＳ３０１）。
【０１０３】
次に、ＣＰＵ  ３０２は、変数xに、適当なm次元の初期ベクトルを記憶させる（ステップＳ３０２）。
【０１０４】
そして、ＣＰＵ  ３０２は、スカラー値A(x)/ρ(x)を計算し、この結果を変数wに記憶させる（ステップＳ３０３）。
【０１０５】
ついで、ＣＰＵ  ３０２は、ベクトル値f(x)を計算し、この結果を変数x'に記憶させる（ステップＳ３０４）。
【０１０６】
さらに、ＣＰＵ  ３０２は、スカラー値A(x)/ρ(x)を計算し、この結果を変数w'に記憶させる（ステップＳ３０５）。
【０１０７】
次に、ＣＰＵ  ３０２は、変数xにx'を記憶させて、変数xを更新する（ステップＳ３０６）。
【０１０８】
そして、ＣＰＵ  ３０２は、変数wにw'を加算して、変数wを更新する（ステップＳ３０７）。
【０１０９】
さらに、ＣＰＵ  ３０２は、変数cに１を加算して、変数cを更新する（ステップＳ３０８）。
【０１１０】
この後、ＣＰＵ  ３０２は、cが所定の繰り返し回数Cを超えたか否かを判断する（ステップＳ３０９）。まだ超えていない場合（ステップＳ３０９；Ｎｏ）、ステップＳ３０４に戻る。
【０１１１】
一方、cが所定の繰り返し回数Cを超えた場合（ステップＳ３０９；Ｙｅｓ）、スカラー値w/(c+1)を計算して、この結果を変数sに記憶させて（ステップＳ３１０）、本処理を終了する。sが積分計算の結果である。
【０１１２】
（ベクトル写像の例）
本実施形態にて用いるベクトル写像fは、m次元ベクトルu = (u1, u2, …, um)と写像f1, f2, …, fmを用いて、
f(u) = (f1(u1), f2(u2), …, fm(um))
のように定義できる。
【０１１３】
図４は、ベクトル写像fを定義する写像fi (1≦i≦m)の例の様子を表すグラフである。以下、本図を参照して説明する。
【０１１４】
各写像fiはカオス写像と呼ばれるものであり、上記のように、正接関数のn (n≧2)倍角の公式から得ることができる。これらの数学的背景や証明については、発明者による以下の論文に開示されている。
Ken Umeno, Superposition of chaotic processes with convergence to Levy's stable law, Physical Review E, vol. 58 No. 2. The American Physical Society, Aug. 1998.
【０１１５】
図中のF2(x)は、パラメータとしてn=2, d=1, a=1を用い、φ6を用いた場合の写像fiの様子である。
【０１１６】
図中のF2*(x)は、パラメータとしてn=2, d=1, a=1を用い、φ5を用いた場合の写像fiの様子である。
【０１１７】
図中のF3(x)は、パラメータとしてn=3, d=1, a=1を用い、φ6を用いた場合の写像fiの様子である。
【０１１８】
図中のF3*(x)は、パラメータとしてn=3, d=1, a=1を用い、φ5を用いた場合の写像fiの様子である。
【０１１９】
これらの関数は、以下のように陽に表現することができる。
F2(x) = 2|x|/(1-x2)
F2*(x) = 1/2(|x|-1/|x|)
F3(x) = x|x2-3|/(1-3x2)
F3*(x) = -x(x2-3)/|(1-3x2)|
【０１２０】
図５は、写像fiを繰り返し適用して得られるスカラー列
t, fi(t), fi(fi(t)), fi(fi(fi(t))), …
の極限分布の密度関数を図示したものである。以下、本図を参照して説明する。
【０１２１】
本図では、d=1と固定し、aを種々変化させた場合の密度関数を描いている。この密度関数は、上記の文献にも示される通り、
ρ(x) = ada|x|a-1/(π(1+d2a|x|2a))
のように、解析的に得ることができる。
【０１２２】
本図からは、いずれのn、φjを採用した場合も、aを変化させると、密度関数ρの形状が激しく変化することがわかる。
【０１２３】
（実験結果１）
図６は、以下のように定義される２種類の１次元被積分関数（スカラー写像）A(・), A'(・)を、無限区間(-∞,+∞)上で、種々のカオス写像f1を用いて数値積分した場合の計算誤差の自乗平均（mean square deviation）を、計算ステップ数N = c+1（cは更新回数）に対してプロットしたグラフである。
A(x) = -4x2 sgn(x)/(π(1+x2)5/2)；
A'(x) = -4|x| sgn(x)/(π(1+x2)5/2)
【０１２４】
これらに対して、上記F3(x)および上記F3*(x)を使用し、各ステップ数Nに対して独立試行を１０００回行って誤差を求め、その平均をとった。
【０１２５】
なお、上記A(・), A'(・)を無限区間(-∞,+∞)上で積分すると、いずれもその結果は０となることが解析的にわかっているため、上記数値積分の結果がそのまま計算誤差に相当することとなる。
【０１２６】
被積分関数A(・)ρ(・)の積分に対してF3(x)を、被積分関数A'(・)ρ(・)の積分に対してF3*(x)を、それぞれ用いた場合の誤差の自乗分散は、ステップ数Nに対して、1/N2のオーダーで急激に減少していく。
【０１２７】
一方、被積分関数A(・)ρ(・)の積分に対してF3*(x)を、被積分関数A'(・)ρ(・)の積分に対してF3(x)を、それぞれ用いた場合の誤差の自乗分散は、ステップ数Nに対して、1/Nのオーダーで、これよりもゆるやかに減少していく。
【０１２８】
したがって、乱数発生に用いるカオス写像を種々用意し、その中から適切なものを選択すれば、通常のモンテカルロ計算に比べて、計算誤差の収束がより高速であることがわかる。
【０１２９】
（実験結果２）
図７は、以下のように定義される２次元関数A(・,・)を、２次元の非有界領域M = (-∞,+∞)×(-∞,+∞)上で、本発明の手法により数値積分した場合の計算誤差の自乗平均を、計算ステップNに対してプロットしたグラフである。
A(x,y)
= ((1-x2)/(1+x2)) ((1-3y2)/(1+y2)3/2) sgn(y)
- sgn(x) (1+x2)1/2 |y| sgn(y) (1+y2)1/2
【０１３０】
写像f1としては、上記のF2(x)とF3(x)を、写像f2としては、上記のF2(x)とF3(x)を、それぞれ使用した４通りの場合について、それぞれ実験結果を図示してある。なお、上記実験と同様、各ステップ数Nに対して独立試行を１０００回行って誤差を求め、その平均をとった。
【０１３１】
積分∫M A(x,y)dxdy = 0が成立することが解析的にわかっているため、上記の場合と同様、数値積分の結果がそのまま計算誤差として得られることになる。
【０１３２】
本図からは、４種類の組み合わせの内、写像f1としてF2(x)を、写像f2としてF3(x)を選択した場合に、誤算の自乗分散はステップ数Nに対してオーダー1/N2で高速に減少していく。これ以外の組み合わせの場合は、オーダー1/Nでゆるやかに減少していく。
【０１３３】
この実験結果からも、乱数発生に用いるカオス写像を種々用意し、その中から適切なものを選択すれば、通常のモンテカルロ計算に比べて、計算誤差の収束がより高速であることがわかる。
【０１３４】
なお、本実施形態での上記のカオス写像fの定義にかえて、
f(u) = (f1(u1, u2, ..., um),
f2(u1, u2, ..., um),
fm(u1, u2, ..., um))；
u = (u1, u2, ..., um)
のような定義を用いることもできる。この場合であっても、fを所定のベクトルに繰り返し適用した結果のベクトル列の極限分布の確率密度が、上記のρの条件を満たすようになっていれば、本発明に適用することができる。
【０１３５】
また、各fiとして、上記のウラム＝フォン・ノイマン写像、キュービック写像、クィンティック写像のほか、チェビシェフ写像や一般化チェビシェフ写像を用いることもできる。
【０１３６】
（計算処理の流れ）
図８は、本発明の計算装置において行われる処理の手順、すなわち、本発明の計算方法の手順の第２の実施形態の流れを示すフローチャートである。以下、本図を参照して説明する。なお、各変数や写像、関数等は、上記実施形態と同様である。なお、本処理は、段落００３４で既に述べたΣを含む式の計算を十分に大きなN（以下では定数Cに対応する。）について行うことに対応するものである。
【０１３７】
まず、ＣＰＵ  ３０２は、変数cに値０を記憶させる（ステップＳ８０１）。
【０１３８】
次に、ＣＰＵ  ３０２は、変数xに、適当なm次元の初期ベクトルを記憶させる（ステップＳ８０２）。
【０１３９】
そして、ＣＰＵ  ３０２は、変数wに０を記憶させる（ステップＳ８０３）。
【０１４０】
さらに、ＣＰＵ  ３０２は、変数xの値を用いてスカラー値A(x)/ρ(x)を計算し、この結果を変数wに加算して、変数wを更新する（ステップＳ８０４）。
【０１４１】
さらに、ＣＰＵ  ３０２は、変数cに１を加算して、変数cを更新する（ステップＳ８０５）。
【０１４２】
この後、ＣＰＵ  ３０２は、cが所定の繰り返し回数Cを超えたか否かを判断する（ステップＳ８０６）。
【０１４３】
まだ超えていない場合（ステップＳ８０６；Ｎｏ）、ＣＰＵ  ３０２は、変数xと、ベクトル写像fとを用いて、m次元ベクトル値f(x)を計算し、この結果を変数xに記憶して変数xを更新して（ステップＳ８０７）、ステップＳ８０４に戻る。
【０１４４】
一方、cが所定の繰り返し回数Cを超えた場合（ステップＳ８０６；Ｙｅｓ）、スカラー値w/cを計算して、この結果を変数sに記憶させて（ステップＳ８０９）、本処理を終了する。sが積分計算の結果である。
【０１４５】
本実施形態では、上記実施形態に比べ、変数を記憶するために必要な領域がより少ない。
【０１４６】
（第３の実施形態）
上記のように、カオス写像fがSuper-efficiencyという性質を満たせば、特に高速にモンテカルロ積分計算が収束する。しかしながら、積分関数Aによっては、Super-efficiencyが満たされるように、ベクトル写像fを選択することが難しい場合がある。本実施形態は、これを自動的に選択するものである。
【０１４７】
図９は、このようなベクトル写像fを選択する選択手順の流れを示すフローチャートである。
【０１４８】
まず、数値積分に用いる際の候補となるベクトル写像を複数用意する（ステップＳ９０１）。たとえば、ベクトル写像fの各成分の計算に用いる写像として、上記のF2、F2*、F3、F3*のように、各種のパラメータを変更したものを用意する。
【０１４９】
次に、候補のベクトル写像のそれぞれについて、上記の多重積分計算方法を、所定の更新回数（繰り返し回数）c1回と、c2 (c1＜c2)回とで、適用し、積分の近似値を得る（ステップＳ９０２）。各候補の収束がどの程度かを見積もるためである。
【０１５０】
次に、得られた複数の数値積分の値を比較して、最も収束の速いベクトル写像を選択する（ステップＳ９０３）。収束の速さは、たとえば、以下のように判定すればよい。
【０１５１】
それぞれのベクトル写像について、c1回の繰り返しにより数値積分の値と、c2回の繰り返しにより数値積分の値との差を求める。次に、この差の絶対値が最も小さいベクトル写像を選択する。これが最も小さい、ということは、それだけ収束が速いことに相当するからである。
【０１５２】
このほか、c2回の計算を省略して、得られた複数の数値積分の値の平均をとり、その平均に最も近い値を出力したベクトル写像を選択してもよい。
【０１５３】
最後に、当該ベクトル写像を用いて、C (c2＜C)回の更新回数で上記の多重積分計算を行い、数値積分の結果を得る（ステップＳ９０４）。
【０１５４】
本実施形態により、収束が高速なSuper-efficiency性を満たすベクトル写像がわからない場合であっても、自動的に収束が高速なベクトル写像を選んで、高速に数値積分を行うことができる。
【０１５５】
（第４の実施形態）
モンテカルロ計算では、乱数の個数（繰り返しの回数、本発明における更新回数）を定めて数値積分を行う場合が多いが、本実施形態では、誤差が所定の値よりも小さくなった場合に、繰り返しを終了する。
【０１５６】
たとえば、上記第２の実施形態の繰り返し計算の各回（１回目、２回目、…）において、スカラー値w/cを計算する。以下、この値を
s1, s2, …, si, …
のように表記する。
【０１５７】
上記実施形態では、繰り返し回数が所定の更新回数Cを超えることをもって、計算の終了としていた。本実施形態では、所定の誤差率εと、最低繰り返し回数L (L≧2)と、を定める。以下の終了条件を満たす場合、
i≧L；
|si−si-1|≦ε|si|
スカラー値siを数値積分の結果として出力するのである。
【０１５８】
図１０は、本実施形態における数値積分の計算方法の手順を示すフローチャートである。各変数や写像、関数は、上記の実施形態と同様であるが、新たに、主記憶装置３０３内に、スカラー変数s'領域を確保する。
【０１５９】
まず、ＣＰＵ  ３０２は、変数cに値０を記憶させる（ステップＳ１００１）。
【０１６０】
次に、ＣＰＵ  ３０２は、変数xに、適当なm次元の初期ベクトルを記憶させる（ステップＳ１００２）。
【０１６１】
そして、ＣＰＵ  ３０２は、変数wに０を記憶させる（ステップＳ１００３）。
【０１６２】
さらに、ＣＰＵ  ３０２は、変数xの値を用いてスカラー値A(x)/ρ(x)を計算し、この結果を変数wに加算して、変数wを更新する（ステップＳ１００４）。
【０１６３】
さらに、ＣＰＵ  ３０２は、変数cに１を加算して、変数cを更新する（ステップＳ１００５）。
【０１６４】
ついで、ＣＰＵ  ３０２は、スカラー値w/cを計算して、この結果を変数s'に記憶させる（ステップＳ１００６）。
【０１６５】
この後、ＣＰＵ  ３０２は、終了条件「c≧Lかつ|s'−s|≦|εs'|」を満たすか否かを判断する（ステップＳ１００７）。
【０１６６】
満たさない場合（ステップＳ１００７；Ｎｏ）、ＣＰＵ  ３０２は、変数sにs'を記憶して更新し（ステップＳ１００８）、変数xと、ベクトル写像fとを用いて、m次元ベクトル値f(x)を計算し、この結果を変数xに記憶して変数xを更新して（ステップＳ１００９）、ステップＳ１００４に戻る。
【０１６７】
一方、終了条件を満たす場合（ステップＳ１００７；Ｙｅｓ）、値s'を変数sに記憶させて（ステップＳ１０１０）、本処理を終了する。sが積分計算の結果である。
【０１６８】
誤差率εは、たとえば０.０１（１パーセント）や０.００１（０.１パーセント）など、所望の任意の値を選択することができる。
【０１６９】
最低繰り返し回数Lは、利用できる計算時間や積分関数Aの種類等に応じて決める。たとえば、１００回〜１０００回程度とすることができる。
【０１７０】
このほか、終了条件が繰り返しの中で、k回（たとえば２〜５０回）連続して満たされた場合、すなわち、
i+k≧L；
|si+1−si|≦ε|si+1|；
|si+2−si+1|≦ε|si+2|；
~~…；
|si+k−si+k-1|≦ε|si+k|
の場合に、繰り返しを終了して、si+kを数値積分の結果として出力する、などの実施態様を採用することもでき、この実施態様も本発明の範囲に含まれる。
【０１７１】
（第５の実施形態）
上記実施例では、非有界領域は(-∞,+∞)×…×(-∞,+∞)であり、全空間についての積分を求めていた。本実施形態は、これに限られない非有界領域Mについて数値積分をする場合の手法である。このような非有界領域Mには、たとえば、(0,+∞)や、(-∞,0)×(0,+∞)や、(0,1)×(-∞,+∞)などが考えられる。
【０１７２】
このようにするためには、第２の実施形態のステップＳ８０７の後に、更新した変数xの値が非有界領域Mに含まれるか否かを判断し、含まれる場合には、ステップＳ８０４に戻り、含まれない場合には、ステップＳ８０７に戻るようにする。
【０１７３】
このようにすることで、任意の非有界領域Mについて、上記の実施形態と同様に、高速に数値積分を行うことができる。
【０１７４】
【発明の効果】
以上説明したように、本発明によれば、多次元の非有界領域上で定義される関数の積分の近似値を高速に計算するのに好適な、多重積分の計算装置、計算方法、および、これらを実現するプログラムを記録した情報記録媒体を提供することができる。
【０１７５】
これらは、金融工学におけるリスク評価、移動体通信や光通信におけるＣＤＭＡの干渉ノイズの評価、インターネットなどの通信におけるトラヒック評価などの数多くの産業分野で利用することができる。
【図面の簡単な説明】
【図１】  本発明の多重積分の計算装置の概要構成を示す模式図である。
【図２】  本発明の多重積分の計算装置が実装される情報処理装置の概要構成を示す模式図である。
【図３】  本発明の多重積分の計算方法の手順を示すフローチャートである。
【図４】  ベクトル写像fを定義する写像fiの例を示すグラフである。
【図５】  写像fiを繰り返し適用して得られるスカラー列の極限分布の密度関数を示すグラフである。
【図６】  １次元被積分関数を無限区間(-∞,+∞)で数値積分した場合の計算ステップ数に対する誤差の自乗平均を示すグラフである。
【図７】  ２次元被積分関数を非有界領域(-∞,+∞)×(-∞,+∞)で数値積分した場合の計算ステップ数に対する誤差の自乗平均を示すグラフである。
【図８】  本発明の多重積分の計算方法の手順を示すフローチャートである。
【図９】  本発明のベクトル写像を選択する手順を示すフローチャートである。
【図１０】  本発明のベクトル写像を選択する手順を示すフローチャートである。
【符号の説明】
１０１  多重積分の計算装置
１１１  第１の記憶部
１１２  第２の記憶部
１２１  第１の計算部
１２２  第２の計算部
１３１  更新部
１４１  出力部
３０１  情報処理装置
３０２  ＣＰＵ
３０３  主記憶装置
３０４  外部記憶装置
３０５  表示装置
３０６  入力装置
３０８  ＲＯＭ
【図１】
0003711223.tif 000002
【図２】
0003711223.tif 000003
【図３】
0003711223.tif 000004
【図４】
0003711223.tif 000005
【図５】
0003711223.tif 000006
【図６】
0003711223.tif 000007
【図７】
0003711223.tif 000008
【図８】
0003711223.tif 000009
【図９】
0003711223.tif 000010
【図１０】
0003711223.tif 000011