図面 (/)

技術 逆ラプラス変換プログラム、逆ラプラス変換のためのテーブル作成プログラム、逆ラプラス変換の数値解算出プログラム、および逆ラプラス変換装置

出願人 国立大学法人京都大学国立大学法人群馬大学
発明者 藤原宏志齋藤三郎松浦勉
出願日 2008年8月6日 (12年4ヶ月経過) 出願番号 2008-203384
公開日 2009年3月26日 (11年9ヶ月経過) 公開番号 2009-064424
状態 拒絶査定
技術分野 複合演算
主要キーワード 正則化法 最良近似 繰上り 逆ラプラス変換 正則化パラメータ 多倍長演算 C言語 観測空間
関連する未来課題
重要な関連分野

この項目の情報は公開日時点(2009年3月26日)のものです。
また、この項目は機械的に抽出しているため、正しく解析できていない場合があります

図面 (20)

課題

原像の導函数が増大する場合、または原像が原点でない場合においても、逆ラプラス変換数値解を計算することができる逆ラプラス変換プログラム、逆ラプラス変換のためのテーブル作成プログラム、逆ラプラス変換の数値解算出プログラム、および逆ラプラス変換装置を提供する。

解決手段

テーブル作成部4は、原点で零であり、かつ絶対連続な函数からなる重み付き再生核ヒルベルト空間上で、重み付き二乗積分空間観測空間としたチホノフ正則化法により導かれる第二種積分方程式離散化して得られる連立方程式の解を求め、連立方程式の解に基づく第二種積分方程式の数値解を含む情報を記述したHテーブルを作成する。逆変換部5は、Hテーブルを参照して、第二種積分方程式の数値解と軟化子函数を乗じたラプラス変換像との重み付き二乗可積分空間での内積数値計算で求める。

概要

背景

逆ラプラス変換は、電気電子工学油田調査画像処理をはじめとする様々な分野で広く応用されている。従来、コンピュータを用いて逆ラプラス変換の数値解を算出するためには、複素積分の数値計算による方法や、ラプラス変換表を用いた方法が用いられていた。しかしながら、このような従来の方法では、丸め誤差または離散化誤差などの計算誤差が生じ、数値的に不安定(ill-conditioned)であるという問題がある。

このような数値的不安定性を解消する手法として、再生核と正則化法に基づいて逆ラプラス変換の数値解を求める方法が提案されている(たとえば、非特許文献1を参照)。この方法では、逆ラプラス変換の数値解が、再生核ヒルベルト空間での最良近似函数として求めることができる。
浦他、"Numerical Real Inversion formulas of the Laplace Transform by Using a Fredholm Integral Equation of the Second Kind"、日本数学会、2007年度年会、応用数学分科会講演アブストラクト、p.69-72

概要

原像の導函数が増大する場合、または原像が原点でない場合においても、逆ラプラス変換の数値解を計算することができる逆ラプラス変換プログラム、逆ラプラス変換のためのテーブル作成プログラム、逆ラプラス変換の数値解算出プログラム、および逆ラプラス変換装置を提供する。テーブル作成部4は、原点で零であり、かつ絶対連続な函数からなる重み付き再生核ヒルベルト空間上で、重み付き二乗積分空間観測空間としたチホノフ正則化法により導かれる第二種積分方程式を離散化して得られる連立方程式の解を求め、連立方程式の解に基づく第二種積分方程式の数値解を含む情報を記述したHテーブルを作成する。逆変換部5は、Hテーブルを参照して、第二種積分方程式の数値解と軟化子函数を乗じたラプラス変換像との重み付き二乗可積分空間での内積を数値計算で求める。

目的

それゆえに、本発明の目的は、原像の導函数が増大する場合、または原像が原点で零でない場合においても、再生核と正則化法を用いて、逆ラプラス変換の数値解を計算することができる逆ラプラス変換プログラム、逆ラプラス変換のためのテーブル作成プログラム、逆ラプラス変換の数値解算出プログラム、および逆ラプラス変換装置を提供することである。

効果

実績

技術文献被引用数
1件
牽制数
0件

この技術が所属する分野

(分野番号表示ON)※整理標準化データをもとに当社作成

ライセンス契約や譲渡などの可能性がある特許掲載中! 開放特許随時追加・更新中 詳しくはこちら

請求項1

与えられた正則化パラメータ軟化子パラメータのもとで、ラプラス変換像に対する原像数値解を算出するために、コンピュータに、原点であり、かつ絶対連続な函数からなる重み付き再生核ヒルベルト空間上で、重み付き二乗積分空間観測空間としたチホノフ正則化法により導かれる第二種積分方程式離散化して得られる連立方程式の解を求め、前記連立方程式の解に基づく前記第二種積分方程式の数値解を含む情報を記述したテーブルを作成するステップと、前記作成したテーブルを記憶部に保存するステップと、前記記憶部のテーブルを参照して、前記第二種積分方程式の数値解と軟化子函数を乗じたラプラス変換像との重み付き二乗可積分空間での内積数値計算で求めることにより、原像の数値解を算出するステップと、前記原像の数値解を出力するステップとを実行させる逆ラプラス変換プログラム

請求項2

正則化パラメータをaとし、軟化子パラメータsとし、任意の関数gのラプラス変換像をLgとし、数値解を算出する原像に対するラプラス変換像をF(p)とし、軟化子函数をRs(p)としたときに、原点で零であり、かつ絶対連続な函数fからなる重み付き再生核ヒルベルト空間Hkのノルムは、式(A1)で表わされ、前記重み付き再生核ヒルベルト空間Hkでの再生核K(t,t′)が式(A2)で表わされ、前記重み付き二乗可積分空間は、L2(R+,u(p)dp)で表わされ、ρ(t)とu(p)の間には、式(A3)の条件が成立し、前記第二種積分方程式は、式(A4)で表わされ、前記内積は、式(A5)で表わされる、請求項1記載の逆ラプラス変換プログラム。

請求項3

ρ(t)は、式(A6)で表わされ、u(p)は、式(A7)で表わされ、Rs(p)は、式(A8)で表わされる、請求項2記載の逆ラプラス変換プログラム。

請求項4

前記コンピュータは、Kビット長(Kは自然数)の汎用レジスタと、直前演算により繰上りまたは繰下がりが生じたか否かを示すフラグを格納するフラグレジスタと、前記フラグを参照して、直前の演算により繰上りまたは繰下がりが生じたか否かに応じた演算を行なう演算器とを備え、前記テーブルを作成するステップおよび前記原像の数値解を算出するステップは、前記数値計算の対象となる変数について、その仮数部をK×Nビット(Nは2以上の自然数)の多倍長で表わし、前記変数の仮数部をKビットごとに区切って、区切られた下位側から順次、前記汎用レジスタを用いて前記演算器に演算を行なわせるステップを含む、請求項3記載の逆ラプラス変換プログラム。

請求項5

前記連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、前記テーブルを作成するステップは、前記行列Aの逆行列A-1を算出して、前記逆行列A-1の要素を前記記憶部に記憶するステップと、計算範囲内のすべてのtについて、前記記憶部に記憶されている前記逆行列A-1の要素の値に基づいて、前記連立方程式の解を算出するステップとを含む、請求項3記載の逆ラプラス変換プログラム。

請求項6

前記連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、前記テーブルを作成するステップは、前記行列Aを上三角行列と下三角行列の積に分解して、前記上三角行列と前記下三角行列の要素を前記記憶部に記憶するステップと、計算範囲内のすべてのtについて、前記記憶部に記憶されている前記上三角行列と前記下三角行列の要素の値に基づいて、前記連立方程式の解を算出するステップとを含む、請求項3記載の逆ラプラス変換プログラム。

請求項7

前記テーブルを作成するステップは、前記正則化パラメータaおよび前記軟化子パラメータsと、その他のパラメータn、U、Lとから、式(A9)〜(A14)に従って、それぞれ変数h、xj、pj、qj、aij、H(pj, t)を計算するステップと、前記連立方程式は、式(A15)で表わされ、式(A16)で表わされる行列Aをコレスキー分解し、式(A15)の解を算出するステップと、式(A17)に従って、式(A15)の解から前記第二種積分方程式の数値解{Hin,a,t:i=0,1,2,・・・,n}を算出し、前記第二種積分方程式の数値解を含む情報を記述したテーブルを作成するステップとを含み、前記原像の数値解を算出するステップは、前記テーブルを参照して、式(A18)に従って、ラプラス変換像F(pj)から、前記原像の数値解fa,s(n)(t)を算出するステップを含む、請求項3記載の逆ラプラス変換プログラム。

請求項8

与えられた正則化パラメータと軟化子パラメータのもとで、ラプラス変換像に対する原像の数値解を算出するためのテーブルを作成するために、コンピュータに、原点で零であり、かつ絶対連続な函数からなる重み付き再生核ヒルベルト空間上で、重み付き二乗可積分空間を観測空間としたチホノフ正則化法により導かれる第二種積分方程式を離散化して得られる連立方程式の解を求め、前記連立方程式の解に基づく前記第二種積分方程式の数値解を含む情報を記述したテーブルを作成するステップと、前記作成したテーブルを記憶部に保存するステップとを実行させるテーブル作成プログラム

請求項9

正則化パラメータをaとし、軟化子パラメータsとし、任意の関数gのラプラス変換像をLgとし、原点で零であり、かつ絶対連続な函数fからなる重み付き再生核ヒルベルト空間Hkのノルムは、式(A19)で表わされ、前記重み付き再生核ヒルベルト空間Hkでの再生核K(t,t′)が式(A20)で表わされ、前記重み付き二乗可積分空間は、L2(R+,u(p)dp)で表わされ、ρ(t)とu(p)の間には、式(A21)の条件が成立し、前記第二種積分方程式は、式(A22)で表わされる、請求項8記載のテーブル作成プログラム。

請求項10

ρ(t)は、式(A23)で表わされ、u(p)は、式(A24)で表わされる、請求項9記載のテーブル作成プログラム。

請求項11

前記コンピュータは、Kビット長(Kは自然数)の汎用レジスタと、直前の演算により繰上りまたは繰下がりが生じたか否かを示すフラグを格納するフラグレジスタと、前記フラグを参照して、直前の演算により繰上りまたは繰下がりが生じたか否かに応じた演算を行なう演算器とを備え、前記テーブルを作成するステップは、前記数値計算の対象となる変数について、その仮数部をK×Nビット(Nは2以上の自然数)の多倍長で表わし、前記変数の仮数部をKビットごとに区切って、区切られた下位側から順次、前記汎用レジスタを用いて前記演算器に演算を行なわせるステップを含む、請求項10記載のテーブル作成プログラム。

請求項12

前記連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、前記テーブルを作成するステップは、前記行列Aの逆行列A-1を算出して、前記逆行列A-1の要素を前記記憶部に記憶するステップと、計算範囲内のすべてのtについて、前記記憶部に記憶されている前記逆行列A-1の要素の値に基づいて、前記連立方程式の解を算出するステップとを含む、請求項10記載のテーブル作成プログラム。

請求項13

前記連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、前記テーブルを作成するステップは、前記行列Aを上三角行列と下三角行列の積に分解して、前記上三角行列と前記下三角行列の要素を前記記憶部に記憶するステップと、計算範囲内のすべてのtについて、前記記憶部に記憶されている前記上三角行列と前記下三角行列の要素の値に基づいて、前記連立方程式の解を算出するステップとを含む、請求項10記載のテーブル作成プログラム。

請求項14

前記テーブルを作成するステップは、前記正則化パラメータaおよび前記軟化子パラメータsと、その他のパラメータn、U、Lとから、式(A25)〜(A30)に従って、それぞれ変数h、xj、pj、qj、aij、H(pj, t)を計算するステップと、前記連立方程式は、式(A31)で表わされ、式(A32)で表わされる行列Aをコレスキー分解し、式(A31)の解を算出するステップと、式(A33)に従って、式(A31)の解から前記第二種積分方程式の数値解{Hin,a,t:i=0,1,2,・・・,n}を算出し、前記第二種積分方程式の数値解を含む情報を記述したテーブルを作成するステップとを含む、請求項10記載のテーブル作成プログラム。

請求項15

請求項8に記載の逆ラプラス変換のためのテーブル作成プログラムによって作成されたテーブルを用いて、ラプラス変換像に対する原像の数値解を作成するために、コンピュータに、請求項8に記載の逆ラプラス変換のためのテーブル作成プログラムで作成されたテーブルが入力され、入力されたテーブルを記憶部に保存するステップと、前記記憶部のテーブルを参照して、前記第二種積分方程式の数値解と軟化子函数を乗じたラプラス変換像との重み付き二乗可積分空間での内積を数値計算で求めることにより、原像の数値解を算出するステップと、前記原像の数値解を出力するステップとを実行させる逆ラプラス変換の数値解算出プログラム

請求項16

請求項9に記載の逆ラプラス変換のためのテーブル作成プログラムによって作成されたテーブルを用いて、ラプラス変換像に対する原像の数値解を作成するために、コンピュータに、請求項9に記載の逆ラプラス変換のためのテーブル作成プログラムで作成されたテーブルが入力され、入力されたテーブルを記憶部に保存するステップと、前記記憶部のテーブルを参照して、前記第二種積分方程式の数値解と軟化子函数Rs(p)を乗じたラプラス変換像F(p)との重み付き二乗可積分空間での内積を数値計算で求めることにより、原像の数値解を算出するステップと、前記内積は、式(A34)で表わされ、前記原像の数値解を出力するステップとを実行させる逆ラプラス変換の数値解算出プログラム。

請求項17

請求項10、12、13に記載の逆ラプラス変換のためのテーブル作成プログラムによって作成されたテーブルを用いて、ラプラス変換像に対する原像の数値解を作成するために、コンピュータに、請求項10、12、13に記載の逆ラプラス変換のためのテーブル作成プログラムで作成されたテーブルが入力され、入力されたテーブルを記憶部に保存するステップと、前記記憶部のテーブルを参照して、前記第二種積分方程式の数値解と軟化子函数Rs(p)を乗じたラプラス変換像F(p)との重み付き二乗可積分空間での内積を数値計算で求めることにより、原像の数値解を算出するステップと、前記軟化子函数Rs(p)は、式(A35)で表わされ、前記内積は、式(A36)で表わされ、前記原像の数値解を出力するステップとを実行させる逆ラプラス変換の数値解算出プログラム。

請求項18

請求項11記載の逆ラプラス変換のためのテーブル作成プログラムによって作成されたテーブルを用いて、逆ラプラス変換の数値解を作成するために、Kビット長(Kは自然数)の汎用レジスタと、直前の演算により繰上りまたは繰下がりが生じたか否かを示すフラグを格納するフラグレジスタと、前記フラグを参照して、直前の演算により繰上りまたは繰下がりが生じたか否かに応じた演算を行なう演算器とを備えたコンピュータに、請求項11記載の逆ラプラス変換のためのテーブル作成プログラムで作成されたテーブルが入力され、入力されたテーブルを記憶部に保存するステップと、前記記憶部のテーブルを参照して、前記第二種積分方程式の数値解と軟化子函数Rs(p)を乗じたラプラス変換像F(p)との重み付き二乗可積分空間での内積を数値計算で求めることにより、原像の数値解を算出するステップと、前記軟化子函数Rs(p)は、式(A37)で表わされ、前記内積は、式(A38)で表わされ、前記原像の数値解を出力するステップとを実行させ、前記原像の数値解を算出するステップは、前記数値計算の対象となる変数について、その仮数部をK×Nビット(Nは2以上の自然数)の多倍長で表わし、前記変数の仮数部をKビットごとに区切って、区切られた下位側から順次、前記汎用レジスタを用いて前記演算器に演算を行なわせるステップを含む、逆ラプラス変換の数値解算出プログラム。

請求項19

請求項14記載の逆ラプラス変換のためのテーブル作成プログラムによって作成されたテーブルを用いて、逆ラプラス変換の数値解を作成するために、コンピュータに、請求項14記載の逆ラプラス変換のためのテーブル作成プログラムで作成されたテーブルが入力され、入力されたテーブルを記憶部に保存するステップと、前記記憶部のテーブルを参照して、式(A39)に従って、前記第二種積分方程式の数値解とラプラス変換像F(pj)から前記原像の数値解fa,s(n)(t)を算出するステップと、前記原像の数値解を出力するステップとを実行させる逆ラプラス変換の数値解算出プログラム。

請求項20

与えられた正則化パラメータと軟化子パラメータのもとで、ラプラス変換像に対する原像の数値解を算出する逆ラプラス変換装置であって、原点で零であり、かつ絶対連続な函数からなる重み付き再生核ヒルベルト空間上で、重み付き二乗可積分空間を観測空間としたチホノフ正則化法により導かれる第二種積分方程式を離散化して得られる連立方程式の解を求め、前記連立方程式の解に基づく前記第二種積分方程式の数値解を含む情報を記述したテーブルを作成するテーブル作成部と、前記作成したテーブルを記憶する記憶部と、前記記憶部のテーブルを参照して、前記第二種積分方程式の数値解と軟化子函数を乗じたラプラス変換像との重み付き二乗可積分空間での内積を数値計算で求めることにより、原像の数値解を算出する逆変換部と、前記原像の数値解を出力する出力部とを備えた逆ラプラス変換装置。

請求項21

正則化パラメータをaとし、軟化子パラメータsとし、任意の関数gのラプラス変換像をLgとし、数値解を算出する原像に対するラプラス変換像をF(p)とし、軟化子函数をRs(p)としたときに、原点で零であり、かつ絶対連続な函数fからなる重み付き再生核ヒルベルト空間Hkのノルムは、式(A40)で表わされ、前記重み付き再生核ヒルベルト空間Hkでの再生核K(t,t′)が式(A41)で表わされ、前記重み付き二乗可積分空間は、L2(R+,u(p)dp)で表わされ、ρ(t)とu(p)の間には、式(A42)の条件が成立し、前記第二種積分方程式は、式(A43)で表わされ、前記内積は、式(A44)で表わされる、請求項20記載の逆ラプラス変換装置。

請求項22

ρ(t)は、式(A45)で表わされ、u(p)は、式(A46)で表わされ、Rs(p)は、式(A47)で表わされる、請求項21記載の逆ラプラス変換装置。

請求項23

前記テーブル作成部および前記逆変換部とは、共通して、Kビット長(Kは自然数)の汎用レジスタと、直前の演算により繰上りまたは繰下がりが生じたか否かを示すフラグを格納するフラグレジスタと、前記フラグを参照して、直前の演算により繰上りまたは繰下がりが生じたか否かに応じた演算を行なう演算器とを備え、前記演算器は、前記数値計算の対象となる変数について、その仮数部をK×Nビット(Nは2以上の自然数)の多倍長で表わし、前記変数の仮数部をKビットごとに区切って、区切られた下位側から順次、前記汎用レジスタを用いて演算を行なう、請求項22記載の逆ラプラス変換装置。

請求項24

前記連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、前記テーブル作成部は、前記行列Aの逆行列A-1を算出して、前記逆行列A-1の要素を前記記憶部に記憶し、計算範囲内のすべてのtについて、前記記憶部に記憶されている前記逆行列A-1の要素の値に基づいて、前記連立方程式の解を算出する、請求項22記載の逆ラプラス変換装置。

請求項25

前記連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、前記テーブル作成部は、前記行列Aを上三角行列と下三角行列の積に分解して、前記上三角行列と前記下三角行列の要素を前記記憶部に記憶し、計算範囲内のすべてのtについて、前記記憶部に記憶されている前記上三角行列と前記下三角行列の要素の値に基づいて、前記連立方程式の解を算出する、請求項22記載の逆ラプラス変換装置。

請求項26

前記テーブル作成部は、前記正則化パラメータaおよび前記軟化子パラメータsと、その他のパラメータn、U、Lとから、式(A48)〜(A53)に従って、それぞれ変数h、xj、pj、qj、aij、H(pj, t)を計算し、前記連立方程式は、式(A54)で表わされ、前記テーブル作成部は、式(A55)で表わされる行列Aをコレスキー分解し、式(A54)の解を算出し、前記テーブル作成部は、さらに式(A56)に従って、式(A54)の解から前記第二種積分方程式の数値解{Hin,a,t:i=0,1,2,・・・,n}を算出し、前記第二種積分方程式の数値解を含む情報を記述したテーブルを作成し、前記逆変換部は、前記テーブルを参照して、式(A57)に従って、ラプラス変換像F(pj)から前記原像の数値解fa,s(n)(t)を算出する、請求項22記載の逆ラプラス変換装置。

技術分野

0001

本発明は、逆ラプラス変換プログラム、逆ラプラス変換のためのテーブル作成プログラム、逆ラプラス変換の数値解算出プログラム、および逆ラプラス変換装置に関し、特に再生核と正則化法用いて逆ラプラス変換の数値解を算出する技術に関する。

背景技術

0002

逆ラプラス変換は、電気電子工学油田調査画像処理をはじめとする様々な分野で広く応用されている。従来、コンピュータを用いて逆ラプラス変換の数値解を算出するためには、複素積分の数値計算による方法や、ラプラス変換表を用いた方法が用いられていた。しかしながら、このような従来の方法では、丸め誤差または離散化誤差などの計算誤差が生じ、数値的に不安定(ill-conditioned)であるという問題がある。

0003

このような数値的不安定性を解消する手法として、再生核と正則化法に基づいて逆ラプラス変換の数値解を求める方法が提案されている(たとえば、非特許文献1を参照)。この方法では、逆ラプラス変換の数値解が、再生核ヒルベルト空間での最良近似函数として求めることができる。
浦他、"Numerical Real Inversion formulas of the Laplace Transform by Using a Fredholm Integral Equation of the Second Kind"、日本数学会、2007年度年会、応用数学分科会講演アブストラクト、p.69-72

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

0004

しかしながら、上記の再生核と正則化法を用いる手法では、原像原点、かつ原像の導函数が増大しない絶対連続函数からなる再生核空間を考えており、原像がこれら2つの条件を満たす場合に対してのみ適用可能である。

0005

それゆえに、本発明の目的は、原像の導函数が増大する場合、または原像が原点で零でない場合においても、再生核と正則化法を用いて、逆ラプラス変換の数値解を計算することができる逆ラプラス変換プログラム、逆ラプラス変換のためのテーブル作成プログラム、逆ラプラス変換の数値解算出プログラム、および逆ラプラス変換装置を提供することである。

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

0006

本発明のある局面に係る逆ラプラス変換プログラムは、与えられた正則化パラメータ軟化子パラメータのもとで、ラプラス変換像に対する原像の数値解を算出するために、コンピュータに、原点で零であり、かつ絶対連続な函数からなる重み付き再生核ヒルベルト空間上で、重み付き二乗積分空間観測空間としたチホノフ正則化法により導かれる第二種積分方程式を離散化して得られる連立方程式の解を求め、連立方程式の解に基づく第二種積分方程式の数値解を含む情報を記述したテーブルを作成するステップと、作成したテーブルを記憶部に保存するステップと、記憶部のテーブルを参照して、第二種積分方程式の数値解と軟化子函数を乗じたラプラス変換像との重み付き二乗可積分空間での内積を数値計算で求めることにより、原像の数値解を算出するステップと、原像の数値解を出力するステップとを実行させる。

0007

好ましくは、正則化パラメータをaとし、軟化子パラメータsとし、任意の関数gのラプラス変換像をLgとし、数値解を算出する原像に対するラプラス変換像をF(p)とし、軟化子函数をRs(p)としたときに、原点で零であり、かつ絶対連続な函数fからなる重み付き再生核ヒルベルト空間Hkのノルムは、式(B1)で表わされ、

0008

0009

重み付き再生核ヒルベルト空間Hkでの再生核K(t,t′)が式(B2)で表わされ、

0010

0011

重み付き二乗可積分空間は、L2(R+,u(p)dp)で表わされ、ρ(t)とu(p)の間には、式(B3)の条件が成立し、

0012

0013

第二種積分方程式は、式(B4)で表わされ、

0014

0015

内積は、式(B5)で表わされる。

0016

0017

好ましくは、ρ(t)は、式(B6)で表わされ、u(p)は、式(B7)で表わされ、Rs(p)は、式(B8)で表わされる。

0018

0019

好ましくは、コンピュータは、Kビット長(Kは自然数)の汎用レジスタと、直前演算により繰上りまたは繰下がりが生じたか否かを示すフラグを格納するフラグレジスタと、フラグを参照して、直前の演算により繰上りまたは繰下がりが生じたか否かに応じた演算を行なう演算器とを備え、テーブルを作成するステップおよび原像の数値解を算出するステップは、数値計算の対象となる変数について、その仮数部をK×Nビット(Nは2以上の自然数)の多倍長で表わし、変数の仮数部をKビットごとに区切って、区切られた下位側から順次、汎用レジスタを用いて演算器に演算を行なわせるステップを含む。

0020

好ましくは、連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、テーブルを作成するステップは、行列Aの逆行列A-1を算出して、逆行列A-1の要素を記憶部に記憶するステップと、計算範囲内のすべてのtについて、記憶部に記憶されている逆行列A-1の要素の値に基づいて、連立方程式の解を算出するステップとを含む。

0021

好ましくは、連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、テーブルを作成するステップは、行列Aを上三角行列と下三角行列の積に分解して、上三角行列と下三角行列の要素を記憶部に記憶するステップと、計算範囲内のすべてのtについて、記憶部に記憶されている上三角行列と下三角行列の要素の値に基づいて、連立方程式の解を算出するステップとを含む。

0022

好ましくは、テーブルを作成するステップは、正則化パラメータaおよび軟化子パラメータsと、その他のパラメータn、U、Lとから、式(B9)〜(B14)に従って、それぞれ変数h、xj、pj、qj、aij、H(pj, t)を計算するステップと、

0023

0024

連立方程式は、式(B15)で表わされ、式(B16)で表わされる行列Aをコレスキー分解し、式(B15)の解を算出するステップと、

0025

0026

式(B17)に従って、式(B15)の解から第二種積分方程式の数値解{Hin,a,t:i=0,1,2,・・・,n}を算出し、第二種積分方程式の数値解を含む情報を記述したテーブルを作成するステップとを含み、

0027

0028

原像の数値解を算出するステップは、テーブルを参照して、式(B18)に従って、ラプラス変換像F(pj)から、原像の数値解fa,s(n)(t)を算出するステップを含む。

0029

0030

本発明のある局面に係るテーブル作成プログラムは、与えられた正則化パラメータと軟化子パラメータのもとで、ラプラス変換像に対する原像の数値解を算出するためのテーブルを作成するために、コンピュータに、原点で零であり、かつ絶対連続な函数からなる重み付き再生核ヒルベルト空間上で、重み付き二乗可積分空間を観測空間としたチホノフ正則化法により導かれる第二種積分方程式を離散化して得られる連立方程式の解を求め、連立方程式の解に基づく第二種積分方程式の数値解を含む情報を記述したテーブルを作成するステップと、作成したテーブルを記憶部に保存するステップとを実行させる。

0031

好ましくは、正則化パラメータをaとし、軟化子パラメータsとし、任意の関数gのラプラス変換像をLgとし、原点で零であり、かつ絶対連続な函数fからなる重み付き再生核ヒルベルト空間Hkのノルムは式(B19)で表わされ、

0032

0033

重み付き再生核ヒルベルト空間Hkでの再生核K(t,t′)が式(B20)で表わされ、

0034

0035

重み付き二乗可積分空間は、L2(R+,u(p)dp)で表わされ、ρ(t)とu(p)の間には、式(B21)の条件が成立し、

0036

0037

第二種積分方程式は、式(B22)で表わされる。

0038

0039

好ましくは、ρ(t)は、式(B23)で表わされ、u(p)は、式(B24)で表わされる。

0040

0041

好ましくは、コンピュータは、Kビット長(Kは自然数)の汎用レジスタと、直前の演算により繰上りまたは繰下がりが生じたか否かを示すフラグを格納するフラグレジスタと、フラグを参照して、直前の演算により繰上りまたは繰下がりが生じたか否かに応じた演算を行なう演算器とを備え、テーブルを作成するステップは、数値計算の対象となる変数について、その仮数部をK×Nビット(Nは2以上の自然数)の多倍長で表わし、変数の仮数部をKビットごとに区切って、区切られた下位側から順次、汎用レジスタを用いて演算器に演算を行なわせるステップを含む。

0042

好ましくは、連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、テーブルを作成するステップは、行列Aの逆行列A-1を算出して、逆行列A-1の要素を記憶部に記憶するステップと、計算範囲内のすべてのtについて、記憶部に記憶されている逆行列A-1の要素の値に基づいて、連立方程式の解を算出するステップとを含む。

0043

好ましくは、連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、テーブルを作成するステップは、行列Aを上三角行列と下三角行列の積に分解して、上三角行列と下三角行列の要素を記憶部に記憶するステップと、計算範囲内のすべてのtについて、記憶部に記憶されている上三角行列と下三角行列の要素の値に基づいて、連立方程式の解を算出するステップとを含む。

0044

好ましくは、テーブルを作成するステップは、正則化パラメータaおよび軟化子パラメータsと、その他のパラメータn、U、Lとから、式(B25)〜(B30)に従って、それぞれ変数h、xj、pj、qj、aij、H(pj, t)を計算するステップと、

0045

0046

連立方程式は、式(B31)で表わされ、式(B32)で表わされる行列Aをコレスキー分解し、式(B31)の解を算出するステップと、

0047

0048

式(B33)に従って、式(B31)の解から第二種積分方程式の数値解{Hin,a,t:i=0,1,2,・・・,n}を算出し、第二種積分方程式の数値解を含む情報を記述したテーブルを作成するステップとを含む。

0049

0050

本発明のある局面に係る逆ラプラス変換の数値解算出プログラムは、上記記載の逆ラプラス変換のためのテーブル作成プログラムによって作成されたテーブルを用いて、ラプラス変換像に対する原像の数値解を作成するために、コンピュータに、上記記載の逆ラプラス変換のためのテーブル作成プログラムで作成されたテーブルが入力され、入力されたテーブルを記憶部に保存するステップと、記憶部のテーブルを参照して、第二種積分方程式の数値解と軟化子函数を乗じたラプラス変換像との重み付き二乗可積分空間での内積を数値計算で求めることにより、原像の数値解を算出するステップと、原像の数値解を出力するステップとを実行させる。

0051

本発明のある局面に係る逆ラプラス変換の数値解算出プログラムは、上記記載の逆ラプラス変換のためのテーブル作成プログラムによって作成されたテーブルを用いて、ラプラス変換像に対する原像の数値解を作成するために、コンピュータに、上記記載の逆ラプラス変換のためのテーブル作成プログラムで作成されたテーブルが入力され、入力されたテーブルを記憶部に保存するステップと、記憶部のテーブルを参照して、第二種積分方程式の数値解と軟化子函数Rs(p)を乗じたラプラス変換像F(p)との重み付き二乗可積分空間での内積を数値計算で求めることにより、原像の数値解を算出するステップと、内積は、式(B34)で表わされ、

0052

0053

原像の数値解を出力するステップとを実行させる。
本発明のある局面に係る逆ラプラス変換の数値解算出プログラムは、上記記載の逆ラプラス変換のためのテーブル作成プログラムによって作成されたテーブルを用いて、ラプラス変換像に対する原像の数値解を作成するために、コンピュータに、上記記載の逆ラプラス変換のためのテーブル作成プログラムで作成されたテーブルが入力され、入力されたテーブルを記憶部に保存するステップと、記憶部のテーブルを参照して、第二種積分方程式の数値解と軟化子函数Rs(p)を乗じたラプラス変換像F(p)との重み付き二乗可積分空間での内積を数値計算で求めることにより、原像の数値解を算出するステップと、軟化子函数Rs(p)は、式(B35)で表わされ、

0054

0055

内積は、式(B36)で表わされ、

0056

0057

原像の数値解を出力するステップとを実行させる。
本発明のある局面に係る逆ラプラス変換の数値解算出プログラムは、上記記載の逆ラプラス変換のためのテーブル作成プログラムによって作成されたテーブルを用いて、逆ラプラス変換の数値解を作成するために、Kビット長(Kは自然数)の汎用レジスタと、直前の演算により繰上りまたは繰下がりが生じたか否かを示すフラグを格納するフラグレジスタと、フラグを参照して、直前の演算により繰上りまたは繰下がりが生じたか否かに応じた演算を行なう演算器とを備えたコンピュータに、上記記載の逆ラプラス変換のためのテーブル作成プログラムで作成されたテーブルが入力され、入力されたテーブルを記憶部に保存するステップと、記憶部のテーブルを参照して、第二種積分方程式の数値解と軟化子函数Rs(p)を乗じたラプラス変換像F(p)との重み付き二乗可積分空間での内積を数値計算で求めることにより、原像の数値解を算出するステップと、軟化子函数Rs(p)は、式(B37)で表わされ、

0058

0059

内積は、式(B38)で表わされ、

0060

0061

原像の数値解を出力するステップとを実行させ、原像の数値解を算出するステップは、数値計算の対象となる変数について、その仮数部をK×Nビット(Nは2以上の自然数)の多倍長で表わし、変数の仮数部をKビットごとに区切って、区切られた下位側から順次、汎用レジスタを用いて演算器に演算を行なわせるステップを含む。

0062

本発明のある局面に係る逆ラプラス変換の数値解算出プログラムは、上記記載の逆ラプラス変換のためのテーブル作成プログラムによって作成されたテーブルを用いて、逆ラプラス変換の数値解を作成するために、コンピュータに、上記記載の逆ラプラス変換のためのテーブル作成プログラムで作成されたテーブルが入力され、入力されたテーブルを記憶部に保存するステップと、記憶部のテーブルを参照して、式(B39)に従って、第二種積分方程式の数値解とラプラス変換像F(pj)から原像の数値解fa,s(n)(t)を算出するステップと、

0063

0064

原像の数値解を出力するステップとを実行させる。
本発明のある局面に係る逆ラプラス変換装置は、与えられた正則化パラメータと軟化子パラメータのもとで、ラプラス変換像に対する原像の数値解を算出する逆ラプラス変換装置であって、原点で零であり、かつ絶対連続な函数からなる重み付き再生核ヒルベルト空間上で、重み付き二乗可積分空間を観測空間としたチホノフ正則化法により導かれる第二種積分方程式を離散化して得られる連立方程式の解を求め、連立方程式の解に基づく第二種積分方程式の数値解を含む情報を記述したテーブルを作成するテーブル作成部と、作成したテーブルを記憶する記憶部と、記憶部のテーブルを参照して、第二種積分方程式の数値解と軟化子函数を乗じたラプラス変換像との重み付き二乗可積分空間での内積を数値計算で求めることにより、原像の数値解を算出する逆変換部と、原像の数値解を出力する出力部とを備える。

0065

好ましくは、正則化パラメータをaとし、軟化子パラメータsとし、任意の関数gのラプラス変換像をLgとし、数値解を算出する原像に対するラプラス変換像をF(p)とし、軟化子函数をRs(p)としたときに、原点で零であり、かつ絶対連続な函数fからなる重み付き再生核ヒルベルト空間Hkのノルムは式(B40)で表わされ、

0066

0067

重み付き再生核ヒルベルト空間Hkでの再生核K(t,t′)が式(B41)で表わされ、

0068

0069

重み付き二乗可積分空間は、L2(R+,u(p)dp)で表わされ、ρ(t)とu(p)の間には、式(B42)の条件が成立し、

0070

0071

第二種積分方程式は、式(B43)で表わされ、

0072

0073

内積は、式(B44)で表わされる。

0074

0075

好ましくは、ρ(t)は、式(B45)で表わされ、u(p)は、式(B46)で表わされ、Rs(p)は、式(B47)で表わされる。

0076

0077

テーブル作成部および逆変換部とは、共通して、Kビット長(Kは自然数)の汎用レジスタと、直前の演算により繰上りまたは繰下がりが生じたか否かを示すフラグを格納するフラグレジスタと、フラグを参照して、直前の演算により繰上りまたは繰下がりが生じたか否かに応じた演算を行なう演算器とを備え、演算器は、数値計算の対象となる変数について、その仮数部をK×Nビット(Nは2以上の自然数)の多倍長で表わし、変数の仮数部をKビットごとに区切って、区切られた下位側から順次、汎用レジスタを用いて演算を行なう。

0078

好ましくは、連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、テーブル作成部は、行列Aの逆行列A-1を算出して、逆行列A-1の要素を記憶部に記憶し、計算範囲内のすべてのtについて、記憶部に記憶されている逆行列A-1の要素の値に基づいて、連立方程式の解を算出する。

0079

好ましくは、連立方程式は、Ax=bの形で表わされ、Aは各要素がtに依存しない係数行列であり、xは連立方程式の解を要素とするベクトルであり、bは各要素がtに依存するベクトルであり、テーブル作成部は、行列Aを上三角行列と下三角行列の積に分解して、上三角行列と下三角行列の要素を記憶部に記憶し、計算範囲内のすべてのtについて、記憶部に記憶されている上三角行列と下三角行列の要素の値に基づいて、連立方程式の解を算出する。

0080

好ましくは、テーブル作成部は、正則化パラメータaおよび軟化子パラメータsと、その他のパラメータn、U、Lとから、式(B48)〜(B53)に従って、それぞれ変数h、xj、pj、qj、aij、H(pj, t)を計算し、

0081

0082

連立方程式は、式(B54)で表わされ、テーブル作成部は、式(B55)で表わされる行列Aをコレスキー分解し、式(B54)の解を算出し、

0083

0084

テーブル作成部は、さらに式(B56)に従って、式(B54)の解から第二種積分方程式の数値解{Hin,a,t:i=0,1,2,・・・,n}を算出し、第二種積分方程式の数値解を含む情報を記述したテーブルを作成し、

0085

0086

逆変換部は、テーブルを参照して、式(B57)に従って、ラプラス変換像F(pj)から原像の数値解fa,s(n)(t)を算出する。

0087

発明の効果

0088

本発明によれば、原像の導函数が増大する場合、または原像が原点で零でない場合においても、再生核と正則化法を用いて、逆ラプラス変換の数値解を計算することができる。

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

0089

以下、本発明に係る実施の形態について図面を参照して説明する。
[第1の実施形態]
まず、第1の実施形態における逆ラプラス変換の基本的なアルゴリズムについて説明する。

0090

(再生核空間理論を用いた逆ラプラス変換)
ある自然な函数空間の函数f(t)のラプラス変換は、式(1)で表わされる。式(1)の逆変換を考えると、通常はBromwich積分による複素解析的な逆変換公式を考えるが、正の実軸上の値のみを用いた逆変換を求めたい場合がある。このような正の実軸上の値のみを用いた逆変換は実逆ラプラス変換と呼ばれる。また、f(t)は原像、F(P)をラプラス変換像と呼ばれる。

0091

0092

式(2)で表わされる正の実軸R+上の有限ノルムを持ち、f(0)=0を満たす絶対連続な函数fからなる重み付き再生核ヒルベルト空間Hkを考える。この重み付き再生核ヒルベルト空間Hkは、原像の導函数f(t)′が増大する場合も包含している。

0093

0094

この重み付き再生核ヒルベルト空間Hkは、式(3)で表わされる再生核K(t,t′)を有する。

0095

0096

このとき、重み付き再生核ヒルベルト空間Hk上の重み付き二乗可積分空間L2(R+,u(p)dp)への線形写像(Lf)(p)pは有界であって、式(4)が成立する。

0097

0098

式(4)から、一般論に従って、以下が成立する。
チホノフ正則化とは、任意のg∈L2(R+,u(p)dp)と任意の正則化パラメータa>0に対して、式(5)の右辺を最小化することである。重み付き再生核ヒルベルト空間Hk上で、重み付き二乗可積分空間L2(R+,u(p)dp)を観測空間としたチホノフ正則化によって、式(6)で表わされる最良近似函数fa(t)が導かれる。

0099

0100

ここで、Ka(・,t)は、式(7)〜(9)で表わされる。

0101

0102

上述の重み付き再生核ヒルベルト空間Hkは、f(0)=0を満たすことを条件としている。f(0)≠0なる函数f(t)についても、上述の重み付き再生核ヒルベルト空間Hkの特性が利用できるように、軟化子函数Rs(p)を用いて、式(6)の最良近似函数fa(t)を修正する。修正最良近似函数fa,s(t)は、式(10)で表わされる。ここで、sは、軟化子パラメータである。

0103

0104

さて、式(9)のtについてラプラス変換をとり、t′を改めてtとおくと、式(11)が得られる。

0105

0106

式(12)および(13)を用いると、式(11)から式(14)が成立する。式(14)は、書き下すと、式(15)のように表わされる。また、式(10)から、式(16)が成立する。

0107

0108

式(14)は、未知の函数が積分の中に現れており、第二種積分方程式と呼ばれる。修正最良近似函数fa,ε(t)は、第二種積分方程式の解と、軟化子函数を乗じたラプラス変換像との重み付き二乗可積分空間での内積である式(16)で表わされる。本発明の実施形態の逆ラプラス変換では、修正最良近似函数fa,ε(t)を原像f(t)の近似解として求める。

0109

(具体的な函数の形)
さて、函数ρ(t)、u(p)、Rs(p)の具体例として、式(17)〜(19)をで表わされるものを用いることができる。

0110

0111

さて、函数ρ(t)が、式(17)で表わされる場合には、式(3)は、式(20)のように表わされる。

0112

0113

さらに、式(20)から、以下の式(21)が成立する。

0114

0115

さらに、式(14)は、式(22)および(23)のように表わされ、式(16)は、式(24)のように表わされる。

0116

0117

(離散化)
式(22)〜(24)をコンピュータで計算するために離散化を行なう。離散化の方法として、ここでは、計算を簡略化するためニストレーム法を用いる。

0118

積分区間の上限U、下限L、積分のサンプル数n、正則化パラメータaを用いて、式(25)〜(28)に従って、以下の変数h、xj、pj、wjを定義する。

0119

0120

式(25)〜式(28)を用いて、式(22)を離散化すると式(29)が得られ、式(24)を離散化すると式(30)が得られる。

0121

0122

したがって、式(29)の連立方程式の解{Hi n,a,t:i=0,1,2,・・・,n}を求めて、それらを式(30)に代入することによって、逆ラプラス変換による原像の数値解fa,s(n)(t)が得られる。

0123

(式(29)の連立方程式の解法
式(29)の連立方程式の解を求める方法として、ここでは、以下の方法を用いる。

0124

まず、式(31)、(32)に示すように以下の変数qj、aijを用いる。このとき、式(33)、(34)が成り立つ。

0125

0126

これらの変数を用いると、式(29)は、式(35)のように変形される。式(35)を行列の形に書き直すと、式(36)のようになり、係数行列Aを式(37)で定義することができる。

0127

0128

係数行列AをLU分解し、前進代入および後退代入することによって、式(36)の解、つまり式(29)の解{Hin,a,t:i=0,1,2,・・・,n}を得ることができる。また、変数qjを用いると、式(30)は、式(38)のように変形することができる。

0129

0130

(構成)
図1は、本発明の実施形態の逆ラプラス変換装置の構成を表わす図である。

0131

図1を参照して、この逆ラプラス変換装置1は、コンピュータを用いて実現されるものであって、入力部2と、メインメモリ6と、CPU3(Central Processing Unit)と、出力部13とを備える。メインメモリ6は、プログラム記憶部7と、変数記憶部8と、Hテーブル記憶部12と、三角行列記憶部9とを含む。CPU3は、テーブル作成部4および逆変換部5として機能する。

0132

入力部2は、サンプル数n、上限値U、下限値L、正則化パラメータa、軟化子パラメータs、ラプラス変換像F(pi)などを入力するものであり、たとえばキーボードなどで実現されている。

0133

出力部13は、逆ラプラス変換による原像の数値解fa,s(n)(t)を外部、たとえば表示装置14に出力する。表示装置14では、逆ラプラス変換による原像の数値解fa,s(n)(t)が文字またはグラフの形で表示される。

0134

プログラム記憶部7は、コンピュータを逆ラプラス変換装置1として機能させるための逆ラプラス変換プログラムが記憶されている。逆ラプラス変換プログラムは、USB(Universal Serial Bus)メモリCD−ROM(Compact Disk - Read Only Memory)などのコンピュータ読取り可能な記録媒体であるリムーマブルメディア99を通じて、外部からインストールすることができる。

0135

テーブル作成部4は、プログラム記憶部7内の逆ラプラス変換プログラムにしたがって、式(22)および式(23)で表わされる第二種積分方程式を離散化して得られる式(35)の連立方程式の解を数値計算で求め、数値計算に基づいて得られる第二種積分方程式の数値解{Hi n,a,t:i=0,1,2,・・・,n}をHテーブルに記述する。テーブル作成部4は、プログラム記憶部7内の逆ラプラス変換プログラムにしたがって、式(37)で示される係数行列Aを上三角行列Uと下三角行列Lの積に分解して、上三角行列Uと下三角行列Lの要素を三角行列記憶部9に記憶させ、計算範囲内のすべてのtについて、三角行列記憶部9に記憶されている上三角行列Uと下三角行列Lの要素の値に基づいて、式(35)の解を算出する。

0136

Hテーブル記憶部12は、テーブル作成部4で作成された図2で示されるようなHテーブルを記憶する。

0137

逆変換部5は、Hテーブルを参照して、式(24)を離散化した式(38)を数値計算することにより、逆ラプラス変換による原像の数値解fa,s(n)(t)を算出する。

0138

変数記憶部8は、テーブル作成部4によるテーブル作成処理および逆変換部5による逆ラプラス変換による原像の数値解を算出するのに必要な数値計算の変数の値を記憶する。

0139

三角行列記憶部9は、テーブル作成部4において、式(35)の連立方程式の解を計算する際に作成された下三角行列Lおよび上三角行列Uを記憶する。

0140

(Hテーブルの作成手順
図3は、第1の実施形態におけるHテーブルを作成する手順を表わすフローチャートである。

0141

まず、外部から、入力部2にサンプル数n、上限値U、下限値L、正則化パラメータaが入力される(ステップS101)。

0142

次に、テーブル作成部4は、h=(U−L)/nを計算して、hの値を変数記憶部8に保存する(ステップS102)。

0143

次に、テーブル作成部4は、j=0、1、2、・・・、nについて、xj=L+j×hを計算して、xjの値を変数記憶部8に保存する(ステップS103)。

0144

次に、テーブル作成部4は、j=0、1、2、・・・、nについて、pj=exp{(π/2)×sinh(xj)}を計算して、pjの値を変数記憶部8に保存する(ステップS104)。

0145

次に、テーブル作成部4は、j=0、1、2、・・・、nについて、qj=pj×cosh(xj)を計算して、qjの値を変数記憶部8に保存する(ステップS105)。

0146

次に、テーブル作成部4は、i=0、1、2、・・・、n、およびj=0、1、2、・・・、nについて、aij=(π/2)×(2/(pi+pj)3)×{1+(pi+pj)+(pi+pj)2/2}×exp(−pj−1/pj)を計算して、aijの値を変数記憶部8に保存する(ステップS106)。

0147

次に、テーブル作成部4は、i=0、1、2、・・・、n、およびj=0、1、2、・・・、nについて係数行列Aの成分を計算する。テーブル作成部4は、i=jのときには、係数行列Aのij成分であるa+qi×aiiを計算し、i≠jのときには、係数行列Aのij成分であるqj×aijを計算する(ステップS107)。

0148

次に、テーブル作成部4は、係数行列AをLU分解して、上三角行列Uおよび下三角行列Lの各要素の値を三角行列記憶部9に保存する(ステップS108)。

0149

次に、テーブル作成部4は、t=t0(初期値)とする(ステップS109)。
次に、テーブル作成部4は、H(pi,t)=(2/pi3)×[1+pi+pi2/2−exp(−tpi)×{1+pi(t+1)+pi2(t+1)2/2}]を計算し、H(pi,t)を変数記憶部8に保存する(ステップS110)。

0150

次に、テーブル作成部4は、三角行列記憶部9に記憶されている上三角行列Uおよび下三角行列Lの要素の値をもとに、前進代入および後退代入によって、式(29)の解{Hin,a,t:i=0,1,2,・・・,n}を求めて、変数記憶部8内に記憶する(ステップS111)。

0151

次に、テーブル作成部4は、変数記憶部8内の式(29)の解{Hin,a,t:i=0,1,2,・・・,n}をHテーブルに記憶する(ステップS112)。

0152

次に、テーブル作成部4は、tがtm(終了値)と等しければ(ステップS113でYES)終了し、tがtm(終了値)と異なるならば(ステップS113でNO)、tをΔtだけ増加させて(ステップS114)、ステップS110からの処理を繰返す。

0153

(Hテーブルを用いた逆ラプラス変換の数値解の算出手順
図4は、第1の実施形態におけるHテーブルを用いて逆ラプラス変換の数値解を算出する手順を表わすフローチャートである。

0154

まず、外部から、入力部2に軟化子パラメータsが入力される(ステップS200)。
次に、逆変換部5は、t=t0に設定する(ステップS201)。

0155

次に、逆変換部5は、j=0、SUM=0に設定する(ステップS202)。
次に、逆変換部5は、HテーブルからHjn,a,tを読出し、変数記憶部8からpj、qj、hを読み出す(ステップS203)。

0156

次に、外部から、入力部2にラプラス変換像F(pj)が入力される(ステップS204)。

0157

次に、逆変換部5は、SUM=SUM+(π/2)×h×F(pj)×(1/s2pj2){exp(−spj)−1}2×pj×Hjn,a,t ×qj×exp(−pj−1/pj)を計算し、計算結果を変数記憶部8に保存する(ステップS205)。

0158

次に、逆変換部5は、jがnと異なるならば(ステップS206でNO)、jを1だけ増加させて(ステップS207)、ステップS203からの処理を繰返す。

0159

逆変換部5は、jがnと等しければ(ステップS206でYES)、逆ラプラス変換値fa,s(n)(t)=SUMとする(ステップS208)。

0160

次に、逆変換部5は、tがtmと等しければ(ステップS209でYES)終了し、tがtm(終了値)と異なるならば(ステップS209でNO)、tをΔtだけ増加させて(ステップS210)、ステップS202からの処理を繰返す。

0161

以上のように、第1の実施形態の逆ラプラス変換装置によれば、原像の導函数が増大する場合、または原像が原点で零でない場合においても、逆ラプラス変換の数値解を計算することができる。また、式(36)の連立方程式の係数行列Aがtに依存しない性質を利用して、係数行列Aを一度だけLU分解し、分解して得られた下三角行列Lと上三角行列Uの各要素の値を記憶しておき、すべてのtについて、この記憶してある下三角行列Lと上三角行列Uの要素の値を用いて、連立方程式の解を算出するので、tが変わるごとに係数行列AをLU分解するのに比べて計算量を削減することができる。

0162

また、連立方程式の解{Hin,a,t:i=0,1,2,・・・,n}は、式(29)を参照すると、ラプラス変換像F(pi)に依存しない。したがって、一度だけ式(36)の解を計算しておき保存しておけば、どのようなF(pi)の逆ラプラス変換を求める場合にも利用することができるという特徴がある。

0163

[第1の実施形態の変形例]
本発明の実施形態では、第二種積分方程式を離散化して得られる連立方程式を構成する係数行列AをLU分解することによって連立方程式を解いたが、これに限定するものではなく、係数行列Aの逆行列A-1を計算して、連立方程式を解いてもよい。

0164

図5は、第1の実施形態の変形例の逆ラプラス変換装置の構成を表わす図である。
図5を参照して、この逆ラプラス変換装置81は、第1の実施形態のラプラス変換装置1と相違する点は、テーブル作成部84と、逆行列記憶部89である。

0165

テーブル作成部84は、プログラム記憶部7内の逆ラプラス変換プログラムにしたがって、式(22)および式(23)で表わされる第二種積分方程式を離散化して得られる式(35)の連立方程式の解を数値計算で求め、数値計算に基づいて得られる第二種積分方程式の数値解{Hi n,a,t:i=0,1,2,・・・,n}をHテーブルに書込む。テーブル作成部84は、プログラム記憶部7内の逆ラプラス変換プログラムにしたがって、式(37)で示される係数行列Aの逆行列A-1を計算して、逆行列A-1の要素の値を逆行列記憶部89に記憶させ、計算範囲内のすべてのtについて、逆行列記憶部89に記憶されている逆行列A-1の要素の値に基づいて、式(35)の解を算出する。

0166

逆行列記憶部89は、テーブル作成部4において、式(35)の連立方程式の解を計算する際に作成された逆行列A-1の各要素の値を記憶する。

0167

(Hテーブルの作成手順)
図6は、第1の実施形態の変形例おけるHテーブルを作成する手順を表わすフローチャートである。

0168

まず、外部から、入力部2にサンプル数n、上限値U、下限値L、正則化パラメータaが入力される(ステップS501)。

0169

次に、テーブル作成部4は、h=(U−L)/nを計算して、hの値を変数記憶部8に保存する(ステップS502)。

0170

次に、テーブル作成部4は、j=0、1、2、・・・、nについて、xj=L+j×hを計算して、xjの値を変数記憶部8に保存する(ステップS503)。

0171

次に、テーブル作成部4は、j=0、1、2、・・・、nについて、pj=exp{(π/2)×sinh(xj)}を計算して、pjの値を変数記憶部8に保存する(ステップS504)。

0172

次に、テーブル作成部4は、j=0、1、2、・・・、nについて、qj=pj×cosh(xj)を計算して、qjの値を変数記憶部8に保存する(ステップS505)。

0173

次に、テーブル作成部4は、i=0、1、2、・・・、n、およびj=0、1、2、・・・、nについて、aij=(π/2)×(2/(pi+pj)3)×{1+(pi+pj)+(pi+pj)2/2}×exp(−pj−1/pj)を計算して、aijの値を変数記憶部8に保存する(ステップS506)。

0174

次に、テーブル作成部4は、i=0、1、2、・・・、n、およびj=0、1、2、・・・、nについて係数行列Aの成分を計算する。テーブル作成部4は、i=jのときには、係数行列Aのij成分であるa+qi×aiiを計算し、i≠jのときには、係数行列Aのij成分であるqj×aijを計算する(ステップS507)。

0175

次に、テーブル作成部4は、係数行列Aの逆行列A-1を計算して、逆行列A-1の各要素の値を逆行列記憶部89に保存する(ステップS508)。

0176

次に、テーブル作成部4は、t=t0(初期値)とする(ステップS509)。
次に、テーブル作成部4は、H(pi,t)=(2/pi3)×[1+pi+pi2/2−exp(−tpi)×{1+pi(t+1)+pi2(t+1)2/2}]を計算し、H(pi,t)を変数記憶部8に保存する(ステップS510)。

0177

次に、テーブル作成部4は、逆行列記憶部89に記憶されている逆行列A-1の要素の値をもとに、式(29)の解{Hin,a,t:i=0,1,2,・・・,n}を求めて、変数記憶部8内に記憶する(ステップS511)。

0178

次に、テーブル作成部4は、変数記憶部8内の式(29)の解{Hin,a,t:i=0,1,2,・・・,n}をHテーブルに記憶する(ステップS512)。

0179

次に、テーブル作成部4は、tがtm(終了値)と等しければ(ステップS513でYES)終了し、tがtm(終了値)と異なるならば(ステップS513でNO)、tをΔtだけ増加させて(ステップS514)、ステップS510からの処理を繰返す。

0180

以上のように、第1の実施形態の変形例の逆ラプラス変換装置によれば、式(36)の連立方程式の係数行列Aがtに依存しない性質を利用して、係数行列Aの逆行列A-1を一度だけ計算して、逆行列A-1の要素の値を記憶しておき、すべてのtについて、この記憶してある逆行列A-1の要素の値を用いて、連立方程式の解を算出するので、tが変わるごとに係数行列Aの逆行列A-1を計算するのに比べて計算量を削減することができる。

0181

[第2の実施形態]
第2の実施形態は、式(29)を変形することによって、係数行列Aの代わりに対称行列A′を用いて、第二種積分方程式を離散化して連立方程式の解を算出する逆ラプラス変換装置に関する。

0182

(式(23)の連立方程式の解法)
第1の実施形態と同様に、式(31)〜(34)で表わされる変数qj、aijを用いる。

0183

第2の実施形態では、式(29)を式(39)のように変形する。

0184

0185

式(40)のように、Hi′n,a,tを定義すると、式(39)は、式(41)のように変形される。式(41)を行列の形に書き直すと、式(42)のようになり、行列A′を式(43)のように定義することができる。

0186

0187

行列A′は、対称行列であり、行列A′をコレスキー分解することによって、式(42)の解を得ることができる。式(42)の解{Hi′n,a,t:i=0,1,2,・・・,n}から、式(44)を用いて、式(29)の解{Hin,a,t:i=0,1,2,・・・,n}が得られる。

0188

0189

(Hテーブルの作成手順)
図7は、第2の実施形態におけるHテーブルを作成する手順を表わすフローチャートである。

0190

まず、外部から、入力部2にサンプル数n、上限値U、下限値L、正則化パラメータaが入力される(ステップS301)。

0191

次に、テーブル作成部4は、h=(U−L)/nを計算して、hの値を変数記憶部8に保存する(ステップS302)。

0192

次に、テーブル作成部4は、j=0、1、2、・・・、nについて、xj=L+j×hを計算して、xjの値を変数記憶部8に保存する(ステップS303)。

0193

次に、テーブル作成部4は、j=0、1、2、・・・、nについて、pj=exp{(π/2)×sinh(xj)}を計算して、pjの値を変数記憶部8に保存する(ステップS304)。

0194

次に、テーブル作成部4は、j=0、1、2、・・・、nについて、qj=pj×cosh(xj)を計算して、qjの値を変数記憶部8に保存する(ステップS305)。

0195

次に、テーブル作成部4は、i=0、1、2、・・・、n、およびj=0、1、2、・・・、nについて、aij=(π/2)×(2/(pi+pj)3)×{1+(pi+pj)+(pi+pj)2/2}×exp(−pj−1/pj)を計算して、aijの値を変数記憶部8に保存する(ステップS306)。

0196

次に、テーブル作成部4は、i=0、1、2、・・・、n、およびj=0、1、2、・・・、nについて係数行列Aの成分を計算する。テーブル作成部4は、i=jのときには、係数行列Aのij成分であるa/qi×aiiを計算し、i≠jのときには、係数行列Aのij成分をaijとする(ステップS307)。

0197

次に、テーブル作成部4は、係数行列Aをコレスキー分解して、下三角行列L、および上三角行列LT(Tは、転置を表わす)の各要素を三角行列記憶部9に保存する(ステップS308)。

0198

次に、テーブル作成部4は、t=t0(初期値)とする(ステップS309)。
次に、テーブル作成部4は、H(pi,t)=(2/pi3)×[1+pi+pi2/2−exp(−tpi)×{1+pi(t+1)+pi2(t+1)2/2}]を計算し、H(pi,t)を変数記憶部8に保存する(ステップS310)。

0199

次に、テーブル作成部4は、三角行列記憶部9に記憶されている下三角行列Lおよび上三角行列LTをもとに、前進代入および後退代入によって、式(41)の解{Hi′n,a,t:i=0,1,2,・・・,n}を求めて、変数記憶部8に記憶する(ステップS311)。

0200

次に、テーブル作成部4は、i=0、1、2、・・・、nについて、Hi n,a,t=Hi′n,a,t/qiを計算して、変数記憶部8に記憶する(ステップS312)。

0201

次に、テーブル作成部4は、変数記憶部8内の式(29)の解である{Hin,a,t:i=0,1,2,・・・,n}をHテーブルに記憶する(ステップS313)。

0202

次に、テーブル作成部4は、tがtm(終了値)と等しければ(ステップS314でYES)終了し、tがtm(終了値)と異なるならば(ステップS314でNO)、tをΔtだけ増加させて(ステップS315)、ステップS310からの処理を繰返す。

0203

以上のように、第2の実施形態の逆ラプラス変換装置によれば、第1の実施形態と同様の効果に加えて、さらに第二種積分方程式を離散化して得られる連立方程式を構成する係数行列Aを対称行列A′に変換し、コレスキー分解を利用して連立方程式の解を算出することができるので、第1の実施形態よりもさらに計算量を削減できる。

0204

[第3の実施形態]
第3の実施形態は、第1の実施形態または第2の実施形態における、連立方程式の解{Hi n,a,t:i=0,1,2,・・・,n}および逆ラプラス変換の数値解fa,s(n)(t)を求める演算を多倍長で行なう装置に関する。以下では、第2の実施形態の演算を多倍長で行なう場合について説明するが、第1の実施形態の演算でも同様である。

0205

(構成)
図8は、第3の実施形態のCPU3の具体的な構成を表わす図である。

0206

図8を参照して、CPU3は、制御部53と、演算器54と、汎用レジスタ群56と、フラグレジスタ57とを備える。図1では、CPU3は、テーブル作成部4および逆変換部5として機能するとして説明した。図1図8の関係を説明すると、テーブル作成部4および逆変換部5の機能を具体的に実現するのが、制御部53と、演算器54と、汎用レジスタ群56と、フラグレジスタ57である。つまり、テーブル作成部4および逆変換部5とは、共通に、制御部53と、演算器54と、汎用レジスタ群56と、フラグレジスタ57とを備える。言い換えると、制御部53と、演算器54と、汎用レジスタ群56と、フラグレジスタ57とは、プログラム記憶部7に格納された逆ラプラス変換プログラムにしたがって、第1または第2の実施形態で説明したHテーブルの作成と、逆ラプラス変換の数値解の算出を行なう。

0207

演算器54は、加算、減算、乗算除算、ビットの論理演算、およびビットシフト演算などを行なう。

0208

汎用レジスタ群56は、16個の汎用レジスタRAX、RBX、RCX、RDX、RSI、RDI、RBP、RSP、R8、R9、R10、R11、R12、R13、R14、R15を含む。図8において、かっこ内は、アセンブラ言語において指定されるレジスタ名前を表わす。各レジスタは、64ビット長である。

0209

図9は、フラグレジスタ57の具体的な構成を表わす図である。
図9を参照して、フラグレジスタ57の最下位ビットに、直前の演算において、汎用レジスタ群56の各汎用レジスタの最上位ビット(64ビット目)において、次の上位ビットへの繰上りまたは繰下がりが生じたか否かを表わすキャリーフラグCFの値が保持される。キャリーフラグCFの値が「1」の場合には、繰上りまたは繰下がりが生じたことを示し、キャリーフラグCFの値が「0」の場合には、繰上りまたは繰下がりが生じなかったことを示す。

0210

演算器54による加算には、繰上がり付き加算と繰上りなし加算とがある。X+Yを繰上り加算する場合には、演算器54は、XとYとキャリーフラグCFの値とを加算し、加算結果に繰上りが生じた場合には、キャリーフラグCFの値を「1」にし、加算結果に繰上りが生じなかった場合には、キャリーフラグCFの値を「0」にする。一方、X+Yを繰上りなし加算する場合には、演算器54は、XとYとを加算し、加算結果に繰上りが生じた場合には、キャリーフラグCFの値を「1」にし、加算結果に繰上りが生じなかった場合には、キャリーフラグCFの値を「0」にする。

0211

演算器54による減算にも、繰下がり付き減算と繰下がりなし減算とがある。X−Yを繰下がり減算する場合には、演算器54は、XからYおよびキャリーフラグCFの値を減算し、減算結果に繰下がりが生じた場合には、キャリーフラグCFの値を「1」にし、減算結果に繰下がりが生じなかった場合には、キャリーフラグCFの値を「0」にする。一方、X−Yを繰下がりなし減算する場合には、演算器54は、XからYを減算し、減算結果に繰下がりが生じた場合には、キャリーフラグCFの値を「1」にし、減算結果に繰下がりが生じなかった場合には、キャリーフラグCFの値を「0」にする。

0212

図8の逆ラプラス変換装置で演算の対象となる数値、たとえば、n、U、L、a、h、s、xj、pj、qj、aij、係数行列Aの要素、係数行列Aをコレスキー分解したLおよびLT、H(pi,t)、Hi′n,a,t、Hin,a,t、F(pj)、S、fa,s(n)(t)などの数値は、浮動小数点形式で表わされる。

0213

図10(a)、(b)、(c)は、図8の逆ラプラス変換装置において演算の対象となる数値の表現を説明するための図である。

0214

図10(a)に示されるように、数値は、浮動小数点形式で表現される。符号部が「0」の場合は正の数を表わし、符号部が「1」の場合は負の数を表わす。仮数部のみが多倍長64nビット(nは2以上の自然数)で表わされる。指数部および符号部は、64ビットで表わされる。

0215

図10(b)、(c)に示されるように、仮数部の整数部の値は、通常は「1」であるが、加算および減算時には、小数点の位置を合わせるために仮数部のビットがシフトされる。仮数部の小数点の位置がKビットだけ左にシフトすると、指数部の値がKだけ減少する。逆に、仮数部の小数点の位置がKビットだけ右にシフトすると、指数部の値がKだけ増加する。

0216

図10(b)、(c)に示されるように、汎用レジスタ群56の64ビットの各汎用レジスタを用いて演算を行なうために、メインメモリ6に要素数nの配列aが確保される。仮数部は、64ビットごとに区切られて、区切られた各64ビットが配列aの各要素に格納される。たとえば、最上位から64個のビットが配列の要素a[0]に格納され、次の64個のビットが配列の要素a[1]に格納される。下位側のビットを格納する配列から、順次演算が行なわれる。

0217

(本発明の実施形態の加算処理
本発明の実施形態における仮数部の加算処理について説明する。

0218

図11は、フラグレジスタ57を備えるCPU3の一例であるAMD64の加算処理の内容を表わす図である。同図は、C言語からadd(c,a,b,n)で呼び出されるアセンブラルーチンを表わす図である。

0219

add(c,a,b,n)は、メインメモリ6内の要素数nの配列aと、メインメモリ6内の要素数nの配列bとを加算して、加算結果をメインメモリ6内の要素数nの配列cに蓄積するルーチンである。

0220

add(c,a,b,n)が呼び出されると、制御部53は、カウンタiの値を保持するためのレジスタ%rcxにnの値を設定し、レジスタ%rsiに配列aへのポインタを設定し、レジスタ%rdxに配列bへのポインタを設定し、レジスタ%rdiに配列cへのポインタを設定する。

0221

行目では、制御部53は、最上位ビットでの繰上りの有無を表わす情報を保持するためのレジスタ%raxの値を「0」(つまり、繰上りなし)に設定するとともに、フラグレジスタ57内のキャリーフラグCFの値を「0」(つまり、繰上りなし)に初期化する。

0222

6行目では、制御部53は、メインメモリ6内のa[i−1]の値をレジスタ%r8にロードする。

0223

7行目では、制御部53は、メインメモリ6内のb[i−1]の値をレジスタ%r10にロードする。

0224

8行目では、演算器54は、繰上り付き加算を実行する。すなわち、演算器54は、レジスタ%r8内のデータとレジスタ%r10内のデータと、フラグレジスタ57内のキャリーフラグCFの値とを加算し、加算結果をレジスタ%r10に保存するとともに、加算により繰上りがあった場合に、フラグレジスタ57内のキャリーフラグCFの値を「1」に設定し、繰上りがなかった場合に、フラグレジスタ57内のキャリーフラグCFの値を「0」に設定する。

0225

9行目では、制御部53は、レジスタ%r10内の加算結果を、メインメモリ6内のc[i−1]の領域に保存する。

0226

10行目では、演算器54は、カウンタiを「1」だけ減らす。
11行目では、制御部53は、カウンタiが「0」の場合には、ループを抜け、カウンタiが「0」以外の場合には、ループの先頭に戻るように制御する。

0227

13行目では、演算器54は、即値「0」と、レジスタ%rax内のデータ(=「0」)と、フラグレジスタ57内のキャリーフラグCFの値とを加算し、加算結果をレジスタ%raxに保存する。したがって、レジスタ%raxには、最上位ビットでの繰上りの有無を表わす情報が保持される。

0228

(参考:キャリーフラグを利用しない加算処理)
次に、比較のために、フラグレジスタ57を備えないCPU3の一例であるAlphaにおける仮数部の加算処理について説明する。

0229

図12は、Alphaの加算処理の内容を表わす図である。同図は、C言語からadd(c,a,b,n)で呼び出されるアセンブラルーチンを表わす図である。

0230

add(c,a,b,n)は、メインメモリ6内の要素数nの配列aと、メインメモリ6内の要素数nの配列bとを加算して、加算結果をメインメモリ6内の要素数nの配列cに蓄積するルーチンである。

0231

add(c,a,b,n)が呼び出されると、制御部53は、カウンタiの値を保持するためのレジスタ$19にnの値を設定し、レジスタ$17に配列aへのポインタを設定し、レジスタ$18に配列bへのポインタを設定し、レジスタ$16に配列cへのポインタを設定する。

0232

2行目では、制御部53は、配列a[n]のアドレスを計算して、それをレジスタ$17に格納する。

0233

3行目では、制御部53は、配列b[n]のアドレスを計算して、それをレジスタ$18に格納する。

0234

4行目では、制御部53は、配列c[n]のアドレスを計算して、それをレジスタ$16に格納する。

0235

6行目では、制御部53は、繰上りの有無を表わす情報を保持するためのレジスタ$0の値を「0」(つまり、繰上りなし)に設定する。

0236

9行目では、制御部53は、メインメモリ6内のa[i−1]の値をレジスタ$1にロードする。

0237

10行目では、制御部53は、メインメモリ6内のb[i−1]の値をレジスタ$2にロードする。

0238

12行目では、演算器54は、加算を実行する。すなわち、演算器54は、レジスタ$1内のデータとレジスタ$2内のデータとを加算し、加算結果をレジスタ$5に保存する。

0239

13行目では、演算器54は、12行目の加算により繰上りがあったか否かを調べる。すなわち、制御部53は、レジスタ$5内のデータ(=加算結果c)がレジスタ$2内のデータ(=b)よりも小さい場合には、繰上りがあったことがわかるので、レジスタ$23に「1」を格納し、それ以外の場合には、繰上りがないことがわかるので、レジスタ$23に「0」を格納する。

0240

15行目では、演算器54は、レジスタ$5内のデータ(=加算結果)と、レジスタ$0内に格納されている1回前のループでの繰上りの有無を表わすデータとを加算して、加算結果をレジスタ$6に格納する。

0241

16行目では、演算器54は、15行目の加算により繰上りがあったか否かを調べる。すなわち、制御部53は、レジスタ$6内のデータがレジスタ$0内のデータよりも小さい場合には、繰上りがあったことがわかるので、レジスタ$0に「1」を格納し、それ以外の場合には、繰上りがないことがわかるので、レジスタ$0に「0」を格納する。

0242

18行目では、制御部53は、レジスタ$6内の加算結果を、メインメモリ6内のc[i−1]の領域に保存する。

0243

19行目では、演算器54は、レジスタ$0内のデータと、レジスタ$23内のデータとの論理和をレジスタ$0に格納することにより、13行目の繰上り判定の結果と16行目の繰上り判定の結果とを統合する。このようにして判定結果が統合できるのは、レジスタ$0内のデータとレジスタ$23内のデータの両方が「1」となる場合は、ありえないからである。

0244

21行目では、制御部53は、配列a[i−2]のアドレスを計算して、それをレジスタ$17に格納する。

0245

22行目では、制御部53は、配列b[i−2]のアドレスを計算して、それをレジスタ$18に格納する。

0246

23行目では、制御部53は、配列c[i−2]のアドレスを計算して、それをレジスタ$16に格納する。

0247

25行目では、演算器54は、カウンタiを「1」だけ減らす。
26行目では、制御部53は、カウンタiが「0」の場合には、ループを抜け、カウンタiが「0」以外の場合には、ループの先頭に戻るように制御する。

0248

(加算処理の比較)
図11図12を比較すると、図11の8行目の繰上り付き加算命令adcqと同等の機能を実現するために、図12では、12行目の加算命令addqと、13行目の比較命令cmplt、15行目の加算命令addqと、16行目の比較命令cmpltと、19行目の論理和命令orの5つの命令が必要となる。以上のことから、本発明の実施形態のようにキャリーフラグCFを利用することにより、それを利用しない場合と比べて、仮数部の加算処理の命令数を削減することができることがわかる。

0249

減算処理
本発明の実施形態における仮数部の減算処理について説明する。

0250

図13は、フラグレジスタ57を備えるCPU3の一例であるAMD64の減算処理の内容を表わす図である。同図は、C言語からsub(c,a,b,n)で呼び出されるアセンブラルーチンを表わす図である。

0251

sub(c,a,b,n)は、メインメモリ6内の要素数nの配列aから、メインメモリ6内の要素数nの配列bを減算して、減算結果をメインメモリ6内の要素数nの配列cに蓄積するルーチンである。

0252

sub(c,a,b,n)が呼び出されると、制御部53は、カウンタiの値を保持するためのレジスタ%rcxにnの値を設定し、レジスタ%rsiに配列aへのポインタを設定し、レジスタ%rdxに配列bへのポインタを設定し、レジスタ%rdiに配列cへのポインタを設定する。

0253

2行目では、制御部53は、最上位ビットでの繰下がりの有無を表わす情報を保持するためのレジスタ%raxの値を「0」(つまり、繰下がりなし)に設定するとともに、フラグレジスタ57内のキャリーフラグCFの値を「0」(つまり、繰下がりなし)に初期化する。

0254

6行目では、制御部53は、メインメモリ6内のa[i−1]の値をレジスタ%r8にロードする。

0255

7行目では、制御部53は、メインメモリ6内のb[i−1]の値をレジスタ%r10にロードする。

0256

8行目では、演算器54は、繰下がり付き減算を実行する。すなわち、演算器54は、レジスタ%r8内のデータから、レジスタ%r10内のデータおよびフラグレジスタ57内のキャリーフラグCFの値とを減算し、減算結果をレジスタ%r10に保存するとともに、減算により繰下がりがあった場合に、フラグレジスタ57内のキャリーフラグCFの値を「1」に設定し、繰下がりがなかった場合に、フラグレジスタ57内のキャリーフラグCFの値を「0」に設定する。

0257

9行目では、制御部53は、レジスタ%r10内の加算結果を、メインメモリ6内のc[i−1]の領域に保存する。

0258

10行目では、演算器54は、カウンタiを「1」だけ減らす。
11行目では、制御部53は、カウンタiが「0」の場合には、ループを抜け、カウンタiが「0」以外の場合には、ループの先頭に戻るように制御する。

0259

13行目では、演算器54は、即値「0」と、レジスタ%rax内のデータ(=「0」)と、フラグレジスタ57内のキャリーフラグCFの値とを加算し、加算結果をレジスタ%raxに保存する。したがって、レジスタ%raxには、最上位ビットでの繰り下がりの有無を表わす情報が保持される。

0260

(参考:キャリーフラグを利用しない減算処理)
次に、比較のために、フラグレジスタ57を備えないCPU3の一例であるAlphaにおける仮数部の減算処理について説明する。

0261

図14は、Alphaの減算処理の内容を表わす図である。同図は、C言語からsub(c,a,b,n)で呼び出されるアセンブラルーチンを表わす図である。

0262

sub(c,a,b,n)は、メインメモリ6内の要素数nの配列aから、メインメモリ6内の要素数nの配列bと減算して、減算結果をメインメモリ6内の要素数nの配列cに蓄積するルーチンである。

0263

sub(c,a,b,n)が呼び出されると、制御部53は、カウンタiの値を保持するためのレジスタ$19にnの値を設定し、レジスタ$17に配列aへのポインタを設定し、レジスタ$18に配列bへのポインタを設定し、レジスタ$16に配列cへのポインタを設定する。

0264

2行目では、制御部53は、配列a[n]のアドレスを計算して、それをレジスタ$17に格納する。

0265

3行目では、制御部53は、配列b[n]のアドレスを計算して、それをレジスタ$18に格納する。

0266

4行目では、制御部53は、配列c[n]のアドレスを計算して、それをレジスタ$16に格納する。

0267

6行目では、制御部53は、繰下がりの有無を表わす情報を保持するためのレジスタ$0の値を「0」(つまり、繰下がりなし)に設定する。

0268

9行目では、制御部53は、メインメモリ6内のa[i−1]の値をレジスタ$1にロードする。

0269

10行目では、制御部53は、メインメモリ6内のb[i−1]の値をレジスタ$2にロードする。

0270

12行目では、演算器54は、減算を実行する。すなわち、演算器54は、レジスタ$1内のデータからレジスタ$2内のデータを減算し、減算結果をレジスタ$5に保存する。

0271

13行目では、演算器54は、12行目の減算により繰下がりがあったか否かを調べる。すなわち、制御部53は、レジスタ$1内のデータ(=a)がレジスタ$2内のデータ(=b)よりも小さい場合には、繰下がりがあったことがわかるので、レジスタ$6に「1」を格納し、それ以外の場合には、繰下がりがないことがわかるので、レジスタ$6に「0」を格納する。

0272

15行目では、演算器54は、レジスタ$5内のデータ(=減算結果)から、レジスタ$0内に格納されている1回前のループでの繰下がりの有無を表わすデータを減算して、減算結果をレジスタ$7に格納する。

0273

16行目では、演算器54は、15行目の減算により繰下がりがあったか否かを調べる。すなわち、制御部53は、レジスタ$5内のデータがレジスタ$0内のデータよりも小さい場合には、繰下がりがあったことがわかるので、レジスタ$0に「1」を格納し、それ以外の場合には、繰下がりがないことがわかるので、レジスタ$0に「0」を格納する。

0274

18行目では、制御部53は、レジスタ$7内の減算結果を、メインメモリ6内のc[i−1]の領域に保存する。

0275

19行目では、演算器54は、レジスタ$0内のデータと、レジスタ$6内のデータとの論理和をレジスタ$0に格納することにより、13行目の繰下がり判定の結果と16行目の繰下がり判定の結果とを統合する。このようにして判定結果が統合できるのは、レジスタ$0内のデータとレジスタ$6内のデータの両方が「1」となる場合は、ありえないからである。

0276

21行目では、制御部53は、配列a[i−2]のアドレスを計算して、それをレジスタ$17に格納する。

0277

22行目では、制御部53は、配列b[i−2]のアドレスを計算して、それをレジスタ$18に格納する。

0278

23行目では、制御部53は、配列c[i−2]のアドレスを計算して、それをレジスタ$16に格納する。

0279

25行目では、演算器54は、カウンタiを「1」だけ減らす。
26行目では、制御部53は、カウンタiが「0」の場合には、ループを抜け、カウンタiが「0」以外の場合には、ループの先頭に戻るように制御する。

0280

(減算処理の比較)
図13図14を比較すると、図13の8行目の繰下がり付き減算命令sbbqと同等の機能を実現するために、図14では、12行目の減算命令subqと、13行目の比較命令cmplt、15行目の減算命令subqと、16行目の比較命令cmpltと、19行目の論理和命令orの5つの命令が必要となる。以上のことから、本発明の実施形態のようにキャリーフラグCFを利用することにより、それを利用しない場合と比べて、加算部の減算処理の命令数を削減することができることがわかる。

0281

乗算処理
まず、本発明の実施形態で用いる乗算の基本的な処理内容を説明する。図15に示すように、要素数4の配列Aと、要素数4の配列Bとの乗算を考える。配列A[0]のデータがa0、配列A[1]のデータがa1、配列A[2]のデータがa2、配列A[3]のデータがa3とする。配列B[0]のデータがb0、配列B[1]のデータがb1、配列B[2]のデータがb2、配列B[3]のデータがb3とする。a0、a1、a2、a3、b0、b1、b2、b3のデータ長は、それぞれ64ビットである。図15に示すように、配列Aと、b0、b1、b2、b3それぞれとの乗算結果は、64×5ビットとなる。配列Aと配列Bとを乗算すると、乗算結果は64×8ビットとなる。

0282

浮動小数点表現の仮数部の積を求めることを考えた場合、このような64×8ビットがすべて必要になるわけではない。下位のビットは、実際上無視しても、計算精度に大きな影響はない。

0283

本発明の実施形態では、要素数nの配列A(つまり64×nビット)と、要素数nの配列B(つまり64×nビット)との乗算では、上位半分の64×nビットのみを求めて、要素数nの配列Cに格納する。

0284

図15で示される配列Aと配列Bの乗算の場合には、上位から64×4ビットだけを求める。この場合、A×b3の計算では、a0×b3の計算は必要だが、a1×b3、a2×b3、a3×b3の計算は省略することができる。一方、A×b0の計算では、a0×b0、a1×b0、a2×b0、a3×b0のすべての計算が必要となる。

0285

図16は、本発明の実施形態における乗算処理の内容を表わす図である。
図16を参照して、この乗算処理は、要素数nの配列A(各要素は64ビット)と、要素数nの配列B(各要素は64ビット)とを乗算して、乗算結果のうち上位半分の64×nビットを、要素数nの配列Cに格納する。

0286

8行目のmul_addは、配列Aの先頭から(n−i)個の要素(A[0]、・・・、A[n−i−1])と、配列B[i]とを乗算して、乗算結果を配列C′[i]から後ろに加えるルーチンである。たとえば、図17に示すように、nが4で、iが3の場合には、B[3](=b3)とA[0](=a0)の乗算だけが行なわれて、乗算結果が配列C′[3]から後ろに加えられる。

0287

(本発明の実施形態のmul_add)
図18は、フラグレジスタ57を備えるCPU3の一例であるAMD64のmul_addの内容を表わす図である。同図は、C言語から、mul_add(c,a,b,n)で呼び出されるアセンブラルーチンを表わす図である。

0288

mul_add(c,a,b,n)は、メインメモリ6内の要素数nの配列aと、メインメモリ6内の変数bとを乗算して、乗算結果をメインメモリ6内の要素数n+1の配列cに加算する命令である。

0289

mul_add(c,a,b,n)が呼び出されると、制御部53は、レジスタ%rdiに配列cへのポインタを設定し、レジスタ%rsiに配列aへのポインタを設定し、レジスタ%rdxに乗数bの値を設定し、カウンタiの値を保持するためのレジスタ%rcxにnの値を設定する。

0290

1行目では、制御部53は、繰上りの有無を表わす情報を保持するためのレジスタ%r9の値を「0」(つまり、繰上りなし)に設定するとともに、フラグレジスタ57内のキャリーフラグCFの値を「0」(つまり、繰上りなし)に初期化する。

0291

2行目では、制御部53は、レジスタ%rdx内の乗数bの値をレジスタ%r8に格納する。

0292

5行目では、制御部53は、メインメモリ6内のa[i−1]の値をレジスタ%raxにロードする。

0293

6行目では、制御部53は、レジスタ%rax内のデータ(=a[i−1])と、レジスタ%r8内のデータ(=b)とを乗算する。乗算結果の上位ビットはレジスタ%rdxに格納され、乗算結果の下位ビットはレジスタ%raxに格納される。

0294

9行目では、演算器54は、繰上りなし加算を実行する。すなわち、演算器54は、レジスタ%r9に格納されている前回のループによる繰上りの有無を表わすデータと、レジスタ%rdxに格納されている乗算結果の上位ビットとを加算して、加算結果をレジスタ%rdxに格納する。この加算では、繰上りが生じないので、フラグレジスタ57内のキャリーフラグCFは「0」となる。

0295

10行目では、演算器54は、繰上りなし加算を実行する。すなわち、演算器54は、レジスタ%raxに格納されている乗算結果の下位ビットと、メインメモリ6内の配列c[i]とを加算して、加算結果をメインメモリ6内の配列c[i]に格納する。この加算によって、繰上りが生じた場合は、キャリーフラグCFは「1」となり、繰上りが生じなかった場合は、キャリーフラグCFは「0」となる。

0296

11行目では、制御部53は、繰上りの有無を表わす情報を保持するためのレジスタ%r9の値を「0」(つまり、繰上りなし)に設定する。

0297

12行目では、演算器54は、繰上り付き加算を実行する。すなわち、演算器54は、フラグレジスタ57内のキャリーフラグCFの値と、レジスタ%rdxに格納されている乗算結果の上位ビットと、メインメモリ6内の配列c[i−1]とを加算して、加算結果をメインメモリ6内の配列c[i−1]に格納するとともに、加算により繰上りがあった場合に、フラグレジスタ57内のキャリーフラグCFの値を「1」に設定し、繰上りがなかった場合に、フラグレジスタ57内のキャリーフラグCFの値を「0」に設定する。

0298

13行目では、演算器54は、繰上り付き加算を実行する。すなわち、演算器54は、即値「0」と、レジスタ%r9内の繰上りの有無を表わす情報と、フラグレジスタ57内のキャリーフラグCFの値とを加算して、加算結果をレジスタ%r9に格納する。この加算では、繰上りが生じないので、フラグレジスタ57内のキャリーフラグCFは「0」となる。

0299

15行目では、演算器54は、カウンタiを「1」だけ減らす。
16行目では、制御部53は、カウンタiが「0」の場合には、ループを抜け、カウンタiが「0」以外の場合には、ループの先頭に戻るように制御する。

0300

(参考:キャリーフラグを利用しないmul_add)
次に、比較のために、フラグレジスタ57を備えないCPU3の一例であるAlphaにおけるmul_addについて説明する。

0301

図19は、Alphaのmul_addの内容を表わす図である。同図は、C言語から、mul_add(c,a,b,n)で呼び出されるアセンブラルーチンを表わす図である。

0302

mul_add(c,a,b,n)は、メインメモリ6内の要素数nの配列aと、メインメモリ6内の変数bとを乗算して、乗算結果をメインメモリ6内の要素数n+1の配列cに加算するルーチンである。

0303

mul_add(c,a,b,n)が呼び出されると、制御部53は、カウンタiの値を保持するためのレジスタ$19にnの値を設定し、レジスタ$17に配列aへのポインタを設定し、レジスタ$18に乗数bの値を設定し、レジスタ$16に配列cへのポインタを設定する。

0304

1行目では、制御部53は、繰上がりの有無を表わす情報を保持するためのレジスタ$23の値を「0」(つまり、繰上がりなし)に設定する。

0305

2行目では、制御部53は、配列c[n]のアドレスを計算して、それをレジスタ$16に格納する。

0306

3行目では、制御部53は、配列a[n]のアドレスを計算して、それをレジスタ$17に格納する。

0307

6行目では、制御部53は、メインメモリ6内のa[i−1]の値をレジスタ$1にロードする。

0308

7行目では、制御部53は、メインメモリ6内のc[i]をレジスタ$21にロードする。

0309

8行目では、制御部53は、メインメモリ6内のc[i−1]をレジスタ$20にロードする。

0310

10行目では、演算器54は、レジスタ$1内のa[i−1]とレジスタ$18内のbとを乗算して、乗算結果の下位64ビットの値をレジスタ$2に格納する。

0311

11行目では、演算器54は、レジスタ$1内のa[i−1]とレジスタ$18内のbとを乗算して、乗算結果の上位64ビットの値をレジスタ$1に格納する。

0312

13行目では、演算器54は、レジスタ$1内の乗算結果の上位64ビットと、レジスタ$23に格納されている1回前のループでの繰上りの有無を表わすデータとを加算して、加算結果をレジスタ$1に格納する。

0313

15行目では、演算器54は、レジスタ$2内の乗算結果の下位64ビットと、レジスタ$21内のc[i]とを加算して、レジスタ$21に格納する。

0314

16行目では、制御部53は、15行目の加算により繰上がりがあったか否かを調べる。すなわち、制御部53は、レジスタ$21内のデータ(=c[i])がレジスタ$2内のデータよりも小さい場合には、繰上がりがあったことがわかるので、レジスタ$2に「1」を格納し、それ以外の場合には、繰上がりがないことがわかるので、レジスタ$2に「0」を格納する。

0315

18行目では、演算器54は、レジスタ$1内の乗算結果の上位64ビット(13行目で繰上り処理されたもの)と、レジスタ$20内のc[i−1]とを加算して、レジスタ$20に格納する。

0316

19行目では、制御部53は、18行目の加算により繰上がりがあったか否かを調べる。すなわち、制御部53は、レジスタ$20内のデータ(=c[i−1])がレジスタ$1内のデータよりも小さい場合には、繰上がりがあったことがわかるので、レジスタ$1に「1」を格納し、それ以外の場合には、繰上がりがないことがわかるので、レジスタ$1に「0」を格納する。

0317

21行目では、演算器54は、レジスタ$2内の15行目の加算の繰上りの有無を表わすデータと、レジスタ$20内のc[i−1](上位64ビットが加算後のもの)とを加算して、レジスタ$20に格納する。

0318

22行目では、制御部53は、21行目の加算により繰上がりがあったか否かを調べる。すなわち、制御部53は、レジスタ$20内のデータ(=c[i−1])がレジスタ$2内のデータよりも小さい場合には、繰上がりがあったことがわかるので、レジスタ$2に「1」を格納し、それ以外の場合には、繰上がりがないことがわかるので、レジスタ$2に「0」を格納する。

0319

24行目では、制御部53は、レジスタ$21内のc[i]をメインメモリ6内の領域に保存する。

0320

25行目では、制御部53は、レジスタ$20内のc[i−1]をメインメモリ6内の領域に保存する。

0321

26行目では、演算器54は、レジスタ$1内のデータと、レジスタ$2内のデータとの論理和をレジスタ$23に格納することにより、19行目の繰上がり判定の結果と22行目の繰上がり判定の結果とを統合する。このようにして判定結果が統合できるのは、レジスタ$1内のデータとレジスタ$2内のデータの両方が「1」となる場合は、ありえないからである。

0322

28行目では、制御部53は、配列a[i−2]のアドレスを計算して、それをレジスタ$17に格納する。

0323

29行目では、制御部53は、配列c[i−2]のアドレスを計算して、それをレジスタ$16に格納する。

0324

31行目では、演算器54は、カウンタiを「1」だけ減らす。
32行目では、制御部53は、カウンタiが「0」の場合には、ループを抜け、カウンタiが「0」以外の場合には、ループの先頭に戻るように制御する。

0325

(mul_addおよび乗算処理の比較)
図18図19とを比較すると、図18では、ループ内部の命令数が7であるのに対して、図19では、ループ内部の命令数が18である。以上のことから、本発明の実施形態のようにキャリーフラグCFを利用することにより、それを利用しない場合と比べて、配列と配列の1要素との乗算であるmul_addの命令数を削減することができることがわかる。その結果、キャリーフラグCFを利用することにより、仮数部を表わす配列と配列の乗算の命令数も削減できる。

0326

除算処理
本発明の実施形態では、古典的なニュートン法により除算を行なう。x=a/bを求めるために、まず、1/bを求め、その結果にaを乗ずる。1/bを求めるために、f(x)=(1/x)−bについてニュートン法を適用する。このとき、ニュートン法の漸化式は、たとえば、以下の式(45)、(46)で与えられる。

0327

x0=1 ・・・(45)
xn+1=xn−f(xn)/f′(xn)=xn×(2−b×xn) ・・・(46)
図20は、ニュートン法を用いた除算処理の内容を説明するための図である。図21は、図20中のサブルーチンcmpの内容を表わす図である。

0328

図20において、21行目は、q=b×xnの乗算を実行するルーチンである。
22行目は、q=(2−q)の減算を実行するルーチンである。

0329

23行目は、q=xn×qの乗算を実行するルーチンである。
25行目および26行目において、xn+1(=q)がxnと等しければ、20行目〜30行目のループを抜ける。

0330

32行目は、xn、すなわち(1/b)にaを乗ずるルーチンである。
(除算処理の比較)
図20の21行目、23行目、32行目の乗算ルーチンmulを実行するCPU3がAMD64の場合では、図16および図18に示されるものであるのに対して、CPU3がAlphaの場合では、図16および図19に示されるものである。両者を比較すると、CPU3がAMD64の方が、命令数が少ないことがわかる。

0331

また、図20の22行目の減算ルーチンsubを実行するCPU3がAMD64の場合では、図13に示されるものであるのに対して、CPU3がAlphaの場合では、図14に示されるものである。両者を比較すると、CPU3がAMD64の方が、命令数が少ないことがわかる。以上のことから、本発明の実施形態のようにキャリーフラグCFを利用することにより、それを利用しない場合と比べて、仮数部の除算処理の命令数を削減することができることがわかる。

0332

(計算の精度についての実験結果)
本発明の実施の形態では、仮数部が64×nビットで表わされるが、nを大きくすることで、fa,s(n)(t)の計算精度に影響する正則化パラメータaの値を小さくすることができる。以下では、正則化パラメータaの大きさによって、fa,s(n)(t)がどのように変化するかの実験結果を説明する。f(t)がステップ函数の場合に、ラプラス変換像F(pi)が解析的に求まるが、これを用いて、実験を行なった。

0333

図22(a)は、正則化パラメータa=10-1、10-4、10-8、10-12のときの、計算結果fa,s(n)(t)を示す図である。図22(b)は、正則化パラメータa=10-100、10-400のときの、計算結果fa,s(n)(t)を示す図である。

0334

これらの図から、正則化パラメータaを小さくするほど、fa,s(n)(t)の計算精度が増加することがわかる。図22内のステップ函数の場合には、正則化パラメータaは、少なくとも10-100ぐらい小さくなくては有効な逆ラプラス変換像が得られないことがわかる。正則化パラメータaの値をこのような小さな値にして、かつ現実的な時間内で逆ラプラス変換像を算出するためには、本発明の実施形態で説明したように、キャリーフラグCFを用いて命令数を大幅に削減した多倍長演算が必要不可欠となる。

0335

また、本発明の実施の形態では、式(6)を修正した式(10)を用いることによって、f(0)が0以外の函数についても、逆ラプラス変換の値を高精度で求めることができるという特徴がある。以下では、軟化子パラメータsの大きさによって、fa,s(n)(t)がどのように変化するかの実験結果を説明する。ここでは、式(47)、(48)で表される函数f(t)を実験に用いた。

0336

0<t<1において、f(t)=1 ・・・(47)
1<tにおいて、f(t)=0 (48)
図23(a)は、s=0.1での計算結果fa,s(n)(t)を示す図であり、図23(b)は、s=0.01での計算結果fa,s(n)(t)を示す図である。

0337

これらの図から、軟化子パラメータsを小さくするほど、fa,s(n)(t)の計算精度が増加することがわかる。軟化子パラメータsを目的に応じて適当な値に設定することによって、f(0)が0以外の函数について、逆ラプラス変換の値を高精度で求めることができる。

0338

[第4の実施形態]
第4の実施形態は、Hテーブルを作成する処理と、Hテーブルを用いて逆ラプラス変換の数値解を算出する処理とを別個の装置で行なう構成に関するものである。

0339

図24は、Hテーブルを作成する逆ラプラス変換用テーブル作成装置の構成を表わす図である。図24の逆ラプラス変換用テーブル作成装置は、コンピュータで実現されるものであり、図24において、図1と同一の構成要素については同一の符号を付している。

0340

図24を参照して、CPU63は、第1〜第3の実施形態と異なり、テーブル作成部4としてのみ機能する。

0341

入力部62は、第1〜第3の実施形態と異なり、Hテーブルの作成に必要なサンプル数n、上限値U、下限値L、正則化パラメータaが入力される。

0342

プログラム記憶部65には、第1〜第3の実施形態と異なり、コンピュータをテーブル作成部4として機能させるための逆ラプラス変換のためのテーブル作成プログラムが記憶されている。逆ラプラス変換のためのテーブル作成プログラムは、USB(Universal Serial Bus)メモリ、CD−ROM(Compact Disk - Read Only Memory)などのコンピュータ読取り可能な記録媒体であるリムーマブルメディア67を通じて、外部からインストールすることができる。

0343

変数記憶部66には、第1〜第3の実施形態と異なり、Hテーブル作成にのみ必要な変数が記憶される。

0344

メインメモリ64内のHテーブル記憶部12に記憶されているHテーブルと、変数記憶部66に記憶されている変数h、pi、qi(i=0,1,2,・・・,n)は、リムーマブルメディア67に出力することができる。これにより、リムーマブルメディア67に記憶されたHテーブルおよび変数を他の装置のメインメモリに送り記憶させることで、他の装置で、Hテーブルおよび変数h、pi、qiを利用することができる。

0345

図25は、逆ラプラス変換の数値解を算出する逆ラプラス変換の数値解算出装置の構成を表わす図である。図25の逆ラプラス変換の数値解算出装置は、コンピュータで実現されるものであり、図25において、図1と同一の構成要素については同一の符号を付している。

0346

図25を参照して、CPU73は、第1〜第3の実施形態と個なり、逆変換部5としてのみ機能する。

0347

入力部72は、第1〜第3の実施形態と異なり、逆ラプラス変換の数値解算出のために必要なラプラス変換像F(pi)と、軟化子パラメータsとが入力される。

0348

プログラム記憶部75には、第1〜第3の実施形態と異なり、コンピュータを逆変換部5および出力部13として機能させるための逆ラプラス変換の数値解算出プログラムが記憶されている。逆ラプラス変換の数値解算出プログラムは、USB(Universal Serial Bus)メモリ、CD−ROM(Compact Disk - Read Only Memory)などのコンピュータ読取り可能な記録媒体であるリムーマブルメディア77を通じて、外部からインストールすることができる。

0349

リムーマブルメディア77に記憶されているHテーブルと、変数h、pi、qi(i=0,1,2,・・・,n)は、それぞれメインメモリ74内のHテーブル記憶部12と、変数記憶部76に出力することができる。これにより、図24の装置で作成されてリムーマブルメディア67に記憶されたHテーブルおよび変数h、pi、qi(i=0,1,2,・・・,n)を、図25の装置で利用することができる。

0350

以上のように、第4の実施形態では、Hテーブルの作成と、逆ラプラス変換の数値解の算出とを別個の装置で行なうことができるので、Hテーブルの作成を高速計算機で行なってユーザに提供し、ユーザは、自身のパソコン上で提供を受けたHテーブルを用いて逆ラプラス変換の数値解を算出を行なうようにすることができる。

0351

(変形例)
本発明は、上記の実施の形態に限定されるものではなく、たとえば以下のような変形例を含む。

0352

(1) テーブルに記述する情報
本発明の第4の実施形態では、図24の逆ラプラス変換用テーブル作成装置が、連立方程式の解{Hi n,a,t:i=0,1,2,・・・,n}を記述したHテーブルと、h、pi、qi(i=0,1,2,・・・,n)とをリムーバブルメディアに出力して、図25の逆ラプラス変換の数値解算出装置でこれを利用したが、これに限定するものではない。たとえば、逆ラプラス変換用テーブル作成装置が、Hテーブルだけをリムーバブルメディアに出力し、逆ラプラス変換の数値解算出装置では、n、U、Lを入力して、h、pi、qjを再計算するものとしてもよい。また、逆ラプラス変換用テーブル作成装置が、h×pi×Hi n,a,t×qiを計算して、計算結果を記述したテーブルとpiとをリムーバブルメディアに出力し、逆ラプラス変換の数値解算出装置で、これらを利用して、式(38)で表わされる計算をより簡単に計算するようにしてもよい。

0353

(2)離散化の手法
本発明の実施形態で説明した式(22)〜(24)を離散化する手法は、様々な方法があり、式(25)〜(36)はその一例にすぎず、その他の方法で行なってもよい。

0354

(3) 64ビットレジスタ
本発明の第3の実施形態では、CPUが64ビット長のレジスタを備える場合を例にして説明したが、これは一例であって、Kビット長(Kは自然数)のレジスタ、たとえば16ビット長、32ビット長、128ビット長などのレジスタであってもよい。

0355

(4)LU分解
本発明の第1の実施形態では、第二種積分方程式を離散化して得られる連立方程式を構成する係数行列AをLU分解したが、これに限定するものではなく、下三角行列と上三角行列の積に分解できるものなら、どのような分解方法であってもよい。

0356

今回開示された実施の形態はすべての点で例示であって制限的なものではないと考えられるべきである。本発明の範囲は上記した説明ではなくて特許請求の範囲によって示され、特許請求の範囲と均等の意味および範囲内でのすべての変更が含まれることが意図される。

図面の簡単な説明

0357

第1の実施形態の逆ラプラス変換装置の構成を表わす図である。
Hテーブルを表わす図である。
第1の実施形態におけるHテーブルを作成する手順を表わすフローチャートである。
第1の実施形態におけるHテーブルを用いて逆ラプラス変換の数値解を算出する手順を表わすフローチャートである。
第1の実施形態の変形例の逆ラプラス変換装置の構成を表わす図である。
第1の実施形態の変形例におけるHテーブルを作成する手順を表わすフローチャートである。
第2の実施形態におけるHテーブルを作成する手順を表わすフローチャートである。
第3の実施形態のCPUの具体的な構成を表わす図である。
フラグレジスタの具体的な構成を表わす図である。
(a)、(b)、(c)は、図8の逆ラプラス変換装置において演算の対象となる数値の表現を説明するための図である。
フラグレジスタを備えるCPUの一例であるAMD64の加算処理の内容を表わす図である。
Alphaの加算処理の内容を表わす図である。
フラグレジスタを備えるCPUの一例であるAMD64の減算処理の内容を表わす図である。
Alphaの減算処理の内容を表わす図である。
本発明の実施形態で用いる乗算処理の内容を説明するための図である。
本発明の実施形態における乗算処理の内容を表わす図である。
本発明の実施形態で用いる乗算処理の内容を説明するための図である。
フラグレジスタを備えるCPUの一例であるAMD64のmul_addの内容を表わす図である。
Alphaのmul_addの内容を表わす図である。
ニュートン法を用いた除算処理の内容を説明するための図である。
図20中のサブルーチンcmpの内容を表わす図である。
(a)は、正則化パラメータa=10-1、10-4、10-8、10-12のときの、計算結果fa,s(n)(t)を示す図である。(b)は、正則化パラメータa=10-100、10-400のときの、計算結果fa,s(n)(t)を示す図である。
(a)は、s=0.1での計算結果fa,s(n)(t)を示す図であり、(b)は、s=0.01での計算結果fa,s(n)(t)を示す図である。
Hテーブルを作成する逆ラプラス変換用テーブル作成装置の構成を表わす図である。
逆ラプラス変換の数値解を算出する逆ラプラス変換の数値解算出装置の構成を表わす図である。

符号の説明

0358

1,81逆ラプラス変換装置、61逆ラプラス変換用テーブル作成装置、71 逆ラプラス変換数値解算出装置、2,62,72 入力部、3,63 CPU、4,84テーブル作成部、5逆変換部、6,64,74メインメモリ、7,65,75プログラム記憶部、8,66,76変数記憶部、9三角行列記憶部、12 Hテーブル記憶部、13 出力部、14表示装置、67,77,99リムーバブルメディア、53 制御部、54演算器、56汎用レジスタ群、57フラグレジスタ、89逆行列記憶部。

ページトップへ

この技術を出願した法人

この技術を発明した人物

ページトップへ

関連する挑戦したい社会課題

該当するデータがありません

関連する公募課題

該当するデータがありません

ページトップへ

技術視点だけで見ていませんか?

この技術の活用可能性がある分野

分野別動向を把握したい方- 事業化視点で見る -

(分野番号表示ON)※整理標準化データをもとに当社作成

ページトップへ

おススメ サービス

おススメ astavisionコンテンツ

新着 最近 公開された関連が強い技術

この 技術と関連性が強い技術

関連性が強い 技術一覧

この 技術と関連性が強い人物

関連性が強い人物一覧

この 技術と関連する社会課題

該当するデータがありません

この 技術と関連する公募課題

該当するデータがありません

astavision 新着記事

サイト情報について

本サービスは、国が公開している情報(公開特許公報、特許整理標準化データ等)を元に構成されています。出典元のデータには一部間違いやノイズがあり、情報の正確さについては保証致しかねます。また一時的に、各データの収録範囲や更新周期によって、一部の情報が正しく表示されないことがございます。当サイトの情報を元にした諸問題、不利益等について当方は何ら責任を負いかねることを予めご承知おきのほど宜しくお願い申し上げます。

主たる情報の出典

特許情報…特許整理標準化データ(XML編)、公開特許公報、特許公報、審決公報、Patent Map Guidance System データ