(57)【要約】      （修正有）【課題】実験値をシミュレーション計算により再現するためのパラメータ決定に要する計算量の削減、および、所望の精度を満たすパラメータ決定に要する計算量の削減をする。【解決手段】高い確度でパラメータを決定できる解析対象物を優先的に選択し、計算順序を制御するリストに格納し、そのリストに基づいて未決定原始パラメータを逐次決定し、すべての原始パラメータが決定するまで分子の選択を行う反復計算を繰り返し実行する。【選択図】図２８
2007323099.tif 000002
【特許請求の範囲】
【請求項１】
複数種類の分子、材料、部品、あるいは構造物のいずれかの解析対象物について測定された特性値や物性値の実験値を、シミュレーション計算により再現するために必要な要素パラメータの決定方法であって、
各解析対象物は複数の一つ又は複数の要素から構成され、各要素は複数種類ある要素タイプのうちの一つに分類され、また、各要素は一つ又は複数の要素属性値を持ち、
(1)入力ステップにおいて、
(1-1)シミュレーションにより得られる計算値と実験値の差に対して許容できる最大値である許容誤差εと、
(1-2)計算値の収束を判定する基準値である収束閾値τと、
(1-3)それぞれの解析対象物の計算に必要な情報、またそれぞれの解析対象物の実験値の実験値と、
(1-4)解析対象物を構成する各要素の要素タイプと、解析対象物を構成する各要素の要素属性値を記憶手段へ入力し、
(2)前処理ステップにおいて、
(2-1)要素タイプごとの要素パラメータ値を初期化し、
(2-2)全ての解析対象物に関する未決定要素パラメータ数を含むデータセットがあり、そのデータセットの並び順を未決定要素パラメータ数の昇順に並べ替え、
(3)計算順序リスト生成ステップにおいて、
(3-1)未決定の要素パラメータを含む解析対象物から、未決定の要素パラメータ数が１以下である解析対象物を選択し、
(3-2)選択された解析対象物の前記(2-2)で並べ替えた並び順を計算順序とする計算順序リストを生成し、
(4)一変数方程式求解ステップにおいて、
(4-1)計算順序リストの最上位の解析対象物が未決定要素パラメータを含む場合にのみ、計算値の変化量が収束閾値未満となるまで未決定パラメータを決定するための一変数方程式の求解を反復し、未決定要素パラメータを決定し、計算順序リストの最上位の解析対象物が未決定要素パラメータを含まない場合には下記(4-2)のステップへ進み、
(4-2)決定した要素パラメータを用いて該解析対象物の計算値と実験値の差である計算誤差を算出し、
(4-3)該計算誤差の絶対値と前記許容誤差を比較し、
(4-4)該計算誤差が許容誤差を超える場合には、仮決定の要素パラメータの１つを未決定の要素パラメータとして新規に設定し、ステップ(4-1)に戻り、該計算誤差が許容誤差を超えない場合には、下記(4-5) のステップへ進み、
(4-5)該計算誤差が許容誤差以下の場合には計算順序リストの最上位のデータセットを削除し、計算順序リストに最上位の解析対象物が存在する場合にはステップ(4-1)に戻り、計算順序リストが空の場合にはステップ(5)へ進み、
(5)要素パラメータ更新ステップにおいて、
(5-1)前記一変数方程式求解ステップで決定した要素パラメータを用いて、要素属性値の関数としての、対応する要素タイプの要素パラメータを更新し、
(5-2)更新された要素タイプの要素パラメータ値を記述する関数を用いて、該要素タイプに分類されている未決定の要素パラメータの一部または全部を仮決定し、
(5-3)未決定の要素パラメータを含む解析対象物の複数のデータセットを、未決定の要素パラメータ数の昇順に並べ替え、
(5-4)前記(5-3)に引き続いて未決定の要素パラメータ数が等しい解析対象物を要素パラメータ数の降順に並べ替え、
(5-5)未決定の要素パラメータを含む解析対象物のデータセットの各パラメータの全てが決定しているのか又は一部でも決定していないのかを決定し、
(6)未決定の解析対象物が存在するかを判定し、未決定の解析対象物が存在する場合は前記(3)のステップへ戻り、
未決定の解析対象物が存在しない場合には、処理を終了することを特徴とするシミュレーション用パラメータ決定方法。
【請求項２】
複数種類の分子、材料、部品、あるいは構造物のいずれかの解析対象物について測定された特性値や物性値の実験値を、シミュレーション計算により再現するために必要な要素パラメータの決定方法であって、
各解析対象物は複数の一つ又は複数の要素から構成され、各要素は複数種類ある要素タイプのうちの一つに分類され、また、各要素は一つ又は複数の要素属性値を持ち、
(1)入力ステップにおいて、
(1-1)計算値の収束を判定する基準値である収束閾値τと、
(1-2)それぞれの解析対象物の実験値の計算に必要なその解析対象物を構成する要素の三次元座標の情報、またそれぞれの解析対象物の実験値の実験値と、
(1-3)解析対象物を構成する各要素の要素タイプと、解析対象物を構成する各要素の要素属性値を記憶手段へ入力し、
(2)前処理ステップにおいて、
(2-1)要素タイプごとの要素パラメータ値を初期化し、
(2-2)全ての解析対象物に関する未決定要素パラメータ数を含むデータセットがあり、そのデータセットの並び順を未決定要素パラメータ数の昇順に並べ替え、
(3)計算順序リスト生成ステップにおいて、
(3-1)未決定の要素パラメータを含む解析対象物から、未決定の要素パラメータ数が１である解析対象物を選択し、
(3-2)選択された解析対象物の前記(2-2)で並べ替えた並び順を計算順序とする計算順序リストを生成し、
(4)一変数方程式求解ステップにおいて、
(4-1)計算順序リストの最上位の解析対象物について、計算値の変化量が収束閾値未満となるまで未決定パラメータを決定するための一変数方程式の求解を反復し、未決定要素パラメータを決定し、
(4-2) 計算順序リストの最上位のデータセットを削除し、計算順序リストに最上位の解析対象物が存在する場合にはステップ(4-1)に戻り、計算順序リストが空の場合にはステップ(5)へ進み、
(5)要素パラメータ更新ステップにおいて、
(5-1)前記一変数方程式求解ステップで決定した要素パラメータを用いて、要素属性値の関数としての、対応する要素タイプの要素パラメータを更新し、
(5-2)未決定の要素パラメータを含む解析対象物の複数のデータセットを、未決定の要素パラメータ数の昇順に並べ替え、
(5-3)前記(5-2)に引き続いて未決定の要素パラメータ数が等しい解析対象物を要素パラメータ数の降順に並べ替え、
(5-4)未決定の要素パラメータを含む解析対象物のデータセットの各パラメータの全てが決定しているのか又は一部でも決定していないのかを決定し、
(6)未決定の解析対象物が存在するかを判定し、未決定の解析対象物が存在する場合は前記(3)のステップへ戻り、
未決定の解析対象物が存在しない場合には、処理を終了することを特徴とするシミュレーション用パラメータ決定方法。
【請求項３】
複数種類の分子について測定された各溶媒和エネルギーの実験値を、シミュレーション計算により再現するために必要な原子パラメータの決定方法であって、
各分子は複数の原子から構成され、各原子は複数種類ある原子タイプのうちの一つに分類され、また、各原子は一つ又は複数の原子属性値を持ち、
(1)入力ステップにおいて、
(1-1)シミュレーション計算で得られる溶媒和エネルギーの計算値と溶媒和エネルギーの実験値の差に対して許容できる最大値である許容エネルギー誤差εと、
(1-2)溶媒和エネルギー計算値の収束を判定する基準値である収束エネルギー閾値τと、
(1-3)それぞれの分子の溶媒和エネルギーの計算に必要なその分子を構成する原子の三次元座標の情報、またそれぞれの分子の溶媒和エネルギーの実験値と、
(1-4)分子を構成する各原子の原子タイプと、分子を構成する各原子の原子属性値を記憶手段へ入力し、
(2)前処理ステップにおいて、
(2-1)原子タイプごとの原子パラメータ値を初期化し、
(2-2)全ての分子に関する未決定原子パラメータ数を含むデータセットがあり、そのデータセットの並び順を未決定原子パラメータ数の昇順に並べ替え、
(3)計算順序リスト生成ステップにおいて、
(3-1)未決定の原子パラメータを含む分子から、未決定の原子パラメータ数が１以下である分子を選択し、
(3-2)選択された分子の前記(2-2)で並べ替えた並び順を計算順序とする計算順序リストを生成し、
(4)一変数方程式求解ステップにおいて、
(4-1)計算順序リストの最上位の分子が未決定原子パラメータを含む場合にのみ、溶媒和エネルギー計算値の変化量が収束エネルギー閾値未満となるまで未決定パラメータを決定するための一変数方程式の求解を反復し、未決定原子パラメータを決定し、計算順序リストの最上位の分子が未決定原子パラメータを含まない場合には下記(4-2)のステップへ進み、
(4-2)決定した原子パラメータを用いて該分子の計算値と実験値の差である計算誤差を算出し、
(4-3)該計算誤差の絶対値と前記許容エネルギー誤差を比較し、
(4-4)該計算誤差が許容エネルギー誤差を超える場合には、仮決定の原子パラメータの１つを未決定の原子パラメータとして新規に設定し、ステップ(4-1)に戻り、該計算誤差が許容エネルギー誤差を超えない場合には、下記(4-5) のステップへ進み、
(4-5)該計算誤差が許容誤差以下の場合には計算順序リストの最上位のデータセットを削除し、計算順序リストに最上位の分子が存在する場合にはステップ(4-1)に戻り、計算順序リストが空の場合にはステップ(5)へ進み、
(5)原子パラメータ更新ステップにおいて、
(5-1)前記一変数方程式求解ステップで決定した原子パラメータを用いて、原子属性値の関数としての、対応する原子タイプの原子パラメータを更新し、
(5-2)更新された原子タイプの原子パラメータ値を記述する関数を用いて、該原子タイプに分類されている未決定の原子パラメータの一部または全部を仮決定し、
(5-3)未決定の原子パラメータを含む分子の複数のデータセットを、未決定の原子パラメータ数の昇順に並べ替え、
(5-4)前記(5-3)に引き続いて未決定の原子パラメータ数が等しい分子を原子パラメータ数の降順に並べ替え、
(5-5)未決定の原子パラメータを含む分子のデータセットの各パラメータの全てが決定しているのか又は一部でも決定していないのかを決定し、
(6)未決定の分子が存在するかを判定し、未決定の分子が存在する場合は前記(3)のステップへ戻り、
未決定の分子が存在しない場合には、処理を終了することを特徴とするシミュレーション用パラメータ決定方法。
【請求項４】
複数種類の分子について測定された各溶媒和エネルギーの実験値を、シミュレーション計算により再現するために必要な原子パラメータの決定方法であって、
各分子は複数の原子から構成され、各原子は複数種類ある原子タイプのうちの一つに分類され、また、各原子は一つ又は複数の原子属性値を持ち、
(1)入力ステップにおいて、
(1-1)溶媒和エネルギー計算値の収束を判定する基準値である収束エネルギー閾値τと、
(1-2)それぞれの分子の溶媒和エネルギーの計算に必要なその分子を構成する原子の三次元座標の情報、またそれぞれの分子の溶媒和エネルギーの実験値と、
(1-3)分子を構成する各原子の原子タイプと、分子を構成する各原子の原子属性値を記憶手段へ入力し、
(2)前処理ステップにおいて、
(2-1)原子タイプごとの原子パラメータ値を初期化し、
(2-2)全ての分子に関する未決定原子パラメータ数を含むデータセットがあり、そのデータセットの並び順を未決定原子パラメータ数の昇順に並べ替え、
(3)計算順序リスト生成ステップにおいて、
(3-1)未決定の原子パラメータを含む分子から、未決定の原子パラメータ数が１である分子を選択し、
(3-2)選択された分子の前記(2-2)で並べ替えた並び順を計算順序とする計算順序リストを生成し、
(4)一変数方程式求解ステップにおいて、
(4-1)計算順序リストの最上位の分子について、溶媒和エネルギー計算値の変化量が収束エネルギー閾値未満となるまで未決定パラメータを決定するための一変数方程式の求解を反復し、未決定原子パラメータを決定し、
(4-2) 計算順序リストの最上位のデータセットを削除し、計算順序リストに最上位の分子が存在する場合にはステップ(4-1)に戻り、計算順序リストが空の場合にはステップ(5)へ進み、
(5)原子パラメータ更新ステップにおいて、
(5-1)前記一変数方程式求解ステップで決定した原子パラメータを用いて、原子属性値の関数としての、対応する原子タイプの原子パラメータを更新し、
(5-2)未決定の原子パラメータを含む分子の複数のデータセットを、未決定の原子パラメータ数の昇順に並べ替え、
(5-3)前記(5-2)に引き続いて未決定の原子パラメータ数が等しい分子を原子パラメータ数の降順に並べ替え、
(5-4)未決定の原子パラメータを含む分子のデータセットの各パラメータの全てが決定しているのか又は一部でも決定していないのかを決定し、
(6)未決定の分子が存在するかを判定し、未決定の分子が存在する場合は前記(3)のステップへ戻り、
未決定の分子が存在しない場合には、処理を終了することを特徴とするシミュレーション用パラメータ決定方法。
【請求項５】
前記計算順序リスト生成ステップにおいて、計算順序を表示し、該計算順序内の一つの分子を選択すると、該選択分子の未決定原子パラメータが分類されている原子タイプに関して、原子電荷と平均結合距離が作る平面における原子の分布図と、該選択分子の原子パラメータの決定状況を分子の三次元構造と分子を構成する各原子上に重ねて描画した球により表現する分子構造図を表示し、これらの表示情報に基づきユーザーが計算順序を変更できる機能を有するユーザーインターフェースを備えることを特徴とする請求項３記載のシミュレーション用パラメータ決定方法。
【請求項６】
前記計算順序リスト生成ステップにおいて、計算順序を表示し、該計算順序内の一つの分子を選択すると、該選択分子の未決定原子パラメータが分類されている原子タイプに関して、原子電荷と平均結合距離が作る平面における原子の分布図と、該選択分子の原子パラメータの決定状況を分子の三次元構造と分子を構成する各原子上に重ねて描画した球により表現する分子構造図を表示し、これらの表示情報に基づきユーザーが計算順序を変更できる機能を有するユーザーインターフェースを備えることを特徴とする請求項４記載のシミュレーション用パラメータ決定方法。
【請求項７】
前記計算順序リスト生成ステップにおいて、計算順序リストが生成された後に、あるいはユーザーが計算順序リスト内の分子を並べ替えあるいは削除を実行した後に、原子パラメータ計算の進行状況の概要として決定、計算中、未決定といった分子の計算状態の集計結果、および計算の経過時間を表示し、また、各分子の原子パラメータ進行状況の概要として各分子の計算状態、原子パラメータ数、未決定原子パラメータ数を表示する機能を有するユーザーインターフェースを備えることを特徴とする請求項３記載のシミュレーション用パラメータ決定方法。
【請求項８】
前記計算順序リスト生成ステップにおいて、計算順序リストが生成された後に、あるいはユーザーが計算順序リスト内の分子を並べ替えあるいは削除を実行した後に、原子パラメータ計算の進行状況の概要として決定、計算中、未決定といった分子の計算状態の集計結果、および計算の経過時間を表示し、また、各分子の原子パラメータ進行状況の概要として各分子の計算状態、原子パラメータ数、未決定原子パラメータ数を表示する機能を有するユーザーインターフェースを備えることを特徴とする請求項４記載のシミュレーション用パラメータ決定方法。
【請求項９】
前記原子パラメータ更新ステップにおいて、各原子タイプの原子パラメータ決定状況として、原子タイプごとに全原子数、原子パラメータが決定した原子数、仮決定の原子数、未決定の原子数を表示し、一つの原子タイプを選択すると、該原子タイプに属する原子の原子パラメータを、原子電荷、平均結合距離、原子パラメータが作る空間内で表示する機能を有するユーザーインターフェースを備えることを特徴とする請求項３記載のシミュレーション用パラメータ決定方法。
【請求項１０】
前記原子パラメータ更新ステップにおいて、各原子タイプの原子パラメータ決定状況として、原子タイプごとに全原子数、原子パラメータが決定した原子数、仮決定の原子数、未決定の原子数を表示し、一つの原子タイプを選択すると、該原子タイプに属する原子の原子パラメータを、原子電荷、平均結合距離、原子パラメータが作る空間内で表示する機能を有するユーザーインターフェースを備えることを特徴とする請求項４記載のシミュレーション用パラメータ決定方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、シミュレーション用パラメータ決定方法に関し、特にシミュレーション計算による実験値の予測に必要なパラメータを決定する方法に関する。
【背景技術】
【０００２】
近年の計算機の発達により、科学分野で計算機シミュレーションはますます多くの機会に利用されている。体内で起こる分子レベルの化学反応のシミュレーションは、その一例であり、特に、タンパク質と低分子化合物の結合反応のシミュレーションは、創薬への応用が期待されている。体内における分子の結合反応では、水が組成の約60%を占めるため、結合に直接関与する分子だけではなく、分子の周囲に存在する水溶媒が結合エネルギーに対して大きく寄与する。そこで、計算機シミュレーションには、水溶媒中での結合エネルギーの高精度計算が求められている。創薬に計算機シミュレーションを用いる場合には、実験値に対する計算値のエネルギー誤差が1.4 kcal/mol以内の計算精度が望まれている。この値は、実験的にタンパク質と低分子化合物の結合を測定した結果得られる結合定数における１０倍の誤差に相当する。
【０００３】
水の効果導入前後のエネルギー差として定義する溶媒和エネルギーの高精度計算には、理論的厳密性が高い計算モデルを使う必要がある。そのようなモデルの一つである第一原理量子化学計算に基づいた連続誘電体モデルと、そこに導入するパラメータを、図１を用いて説明する。このモデルでは、計算対象分子101の原子核と電子をあらわに扱うのに対して、水溶媒102は計算対象分子101を包含する誘電体として扱う。また、図中白抜きで示すように、分子が占有する空間103と、それ以外、つまり、誘電体が占有する空間の境界面104を定義する必要がある。境界面は、分子を構成する全ての原子上に配置された球105の重なりが形成する表面として定義できる。近似的に水素原子上には球を配置しない選択をしてもよい。溶媒和エネルギーは、分子と水溶媒の相互作用の結果、境界面上に存在する電荷に基づいて計算する。そのため、境界面の定義は計算精度に直接影響する。そこで、原子上に配置する球の半径106を溶質と溶媒の境界面を規定するパラメータと考え、この値を決定すれば、溶媒和エネルギーの高精度計算が可能となる。今後、この半径を「原子パラメータ」と呼ぶ。
【０００４】
原子パラメータの決定方法の一つに勾配行列を用いる方法がある。以下に、図２および図３を用いて、勾配行列を用いた原子パラメータ決定手順を説明する。
【０００５】
図２は、勾配行列を用いた原子パラメータ決定方法のフローチャートである。まず、入力装置201から収束エネルギー閾値τ2011、分子の三次元座標情報および溶媒和エネルギーの実験値2012、分子を構成する各原子の化学的類似性に基づく分類である原子タイプ2013を入力する。次に、原子パラメータ初期化ステップ202では、原子タイプごとに設定した原子パラメータの値を初期化する。次に、行列生成ステップ203では、全分子の溶媒和エネルギー、および全分子の、各原子パラメータの変化に対する溶媒和エネルギーの勾配を格納する行列を生成する。次に、行列方程式求解ステップ204では、生成した行列をA,  原子パラメータの更新値ベクトルをx、溶媒和エネルギー実験値ベクトルをbとしたとき、目的関数（数１）
【０００６】
【数１】
2007323099.tif 000003
を最小化するような更新値ベクトルxを求解する。次に、原子パラメータ更新ステップ205では、更新値ベクトルに基づき原子パラメータ値を更新する。次に、原子パラメータ計算終了判定ステップ206では、溶媒和エネルギーの残差平方和の変化量がτ未満であるか判定する。変化量がτ以上の場合は、行列生成ステップ203、行列方程式求解ステップ204、原子パラメータ更新ステップ205を反復する。原子パラメータの変化量がτ未満の場合は、出力装置207より、収束した原子パラメータ値2071を出力する。
【０００７】
図３は、勾配行列を用いた原子パラメータ決定方法における、行列生成ステップ、行列方程式求解ステップ、および原子パラメータ更新ステップにおける処理の詳細である。行列生成ステップ301で実行する処理302は、全分子の溶媒和エネルギー3021、および各原子パラメータの変化に対する全解析対象物の溶媒和エネルギー勾配3022の計算である。計算結果は行列3023に格納する。行列内で○印で示す行列要素3024の一つ一つに要する計算量が、連続誘電体モデルを用いた溶媒和エネルギーの評価一回分に相当する。行列方程式求解ステップ303で実行する処理304は、行列生成ステップ301で生成された行数が分子数、列数が（原子タイプ数＋１）からなる行列3041と、（原子タイプ数＋１）個の要素数を持つ解ベクトル3042と、分子数と同一の要素数を持つ溶媒和エネルギー実験値ベクトル3043で構成される行列方程式から定義できる目的関数gの最小化計算である。行列3041と、溶媒和エネルギー実験値ベクトル3043は既知であり、解ベクトル3042が求解対象の未知量である。原子パラメータ更新ステップ305で実行する処理306は、新しい原子パラメータベクトル3061の、更新前の原子パラメータベクトル3062と行列方程式求解ステップ303で得られた解ベクトルから先頭要素のみを削除したベクトル3063の和による更新である。ベクトル3063の各要素は、原子パラメータ値の次回ループ時の変化量を表している。原子パラメータを更新すると、溶媒和エネルギーや、溶媒和エネルギーの原子パラメータ変化に対する勾配が変化するので、g＜τを満たすまで、行列生成ステップから原子パラメータ更新までのステップまでを反復計算する。
【０００８】
このような原子パラメータの決定方法には二つの問題点がある。
第一の問題点は、行列生成ステップ301において全分子の溶媒和エネルギー計算3021、各原子パラメータの変化に対する全分子の溶媒和エネルギー勾配計算3022に要する計算量は膨大なことである。なぜならば、第一原理量子化学計算では、電子間の相互作用を表すために、原子数のm乗のオーダーで(m=3,4,5,..  計算手法により異なる)計算量が増加する分子積分を計算するからである。一回の溶媒和エネルギー計算に要する計算量を１とすると、一回の溶媒和エネルギー勾配計算に要する計算量も有限差分法では１となる。したがって、行列生成に要する総計算量は、行列3023内に○印で示す行列要素3024の数と等しい。行列生成ステップは、溶媒和エネルギーが収束するまで繰り返すので、パラメータ決定に要する総計算量は、（行列要素数×反復回数）となる。
【０００９】
より一般的に、原子パラメータ決定に用いた分子数をN、分子nに含まれる原子数をA(n)、分子nの溶媒和エネルギー計算に要する計算量をf(A(n))とする。先に述べたように、f(A(n))はA(n)のべき乗関数である。また、分子nに含まれる原子タイプ数をp(n)、計算量をTとする。
【００１０】
勾配行列を用いる従来法では、原子パラメータ変化量が収束エネルギー閾値未満となるまで計算を反復する。反復回数をIとすると、計算量は（数２）
【００１１】
【数２】
2007323099.tif 000004
となる。数２の１行目右辺における和記号以降が一回の行列に生成に要する計算量であり、和記号内部の”1”が溶媒和エネルギー計算に相当し、p(n)が原子パラメータ値の変化に対する溶媒和エネルギーの勾配計算に相当する。
【００１２】
第一の問題点を、具体例を用いて補足説明する。図４は、勾配行列を用いた原子パラメータ決定方法の１０分子への適用例である。図中には、原子タイプ401、得られた原子パラメータ値402、および各分子が含有する原子タイプ403を示す。原子タイプ401としてCH3(C), CH3(N), CH2(CC), CH2(CN),NH2(C),NH2(N)の６種類を定義した。括弧の中に示した原子は、結合相手原子を示す。例えば、原子タイプCH2(CN)は、１つの炭素原子と、１つの窒素原子との結合を有するCH2の炭素原子を意味する。各分子に対応する原子タイプ403は○印で示した。例えば、６番目の分子であるCH3-CH3-NH2では、CH3はCH2に結合しているためCH3(C)の原子タイプに属し、CH2はCH3とNH2に結合しているためCH2(CN)の原子タイプに属し、NH2はCH2に結合しているためNH2(C)の原子タイプに属する。その結果、該分子には３箇所の○印が付与されている。本原子パラメータ決定例における例外として、９番目の分子である(CH3)4-Cの非メチル炭素原子には1.7の固定値を用いた。この原子は、４つのメチル基に周囲を包囲されているため、水分子との接触が困難であり、溶媒和エネルギーへの寄与は軽微と考えられるからである。
【００１３】
原子パラメータの初期値は、CH3(C), CH3(N), CH2(CC), CH2(CN)の四種については1.7、NH2(C), NH2(N)の二種については1.4とした。求解の具体的な方法として、毎回の反復ごとに、以下に示す最急勾配法により解ベクトルを近似的に決定し、原子パラメータを更新した。また、原子パラメータの収束を加速するため、４回の反復に一度は履歴解を用いた原子パラメータを更新した。以下に、その詳細を示す。
１） 毎回の反復ごとの更新方法
イ．各分子nを構成する、対称性の観点から等価でない各原子aについて、溶媒和エネルギーEnのテーラー展開を一次項で打ち切り得られた表式とエネルギー残差znに基づいて、原子パラメータrn,a,tおよび、その重みwn,a,tを計算。（数３）、（数４）。
【００１４】
【数３】
2007323099.tif 000005
【００１５】
【数４】
2007323099.tif 000006
【００１６】
ロ．原子パラメータの更新値rtを各分子nについて得られたrn,a,tの重み付き平均として算出。（数５）
【００１７】
【数５】
2007323099.tif 000007
２） ４回の反復に一度の更新方法
イ．各分子nを構成する、対称性の観点から等価でない各原子aについて、反復回数をI、iを整数としたとき、I=4i-3, 4i-2, 4i-1, 4i時の原子パラメータrn,a,t(I)および残差zn,(I)の履歴を近似する直線zn=Arn,a,t+Bを求め、残差が0となる原子パラメータ値をrn,a,t=-B/Aとして算出。
ロ．各原子タイプtについて、原子パラメータの更新値rtを得られたrn,a,tの平均として算出。分母のデルタ関数に関する和は、分子nを構成する対称性の観点から等価でない原子のうち、原子タイプtに属する原子の数である。（数６）
【００１８】
【数６】
2007323099.tif 000008
【００１９】
図５に、上記に説明した勾配行列を用いた原子パラメータ決定方法の適用により得られた計算量および計算誤差を示す。図中横軸は一台のコンピュータを用いた原子パラメータ決定に要した計算時間501と定義する計算量、縦軸は各分子の絶対計算誤差の平均値と定義する計算誤差502である。ひし形プロット503は計算量と計算誤差の関係を示す。この方法では、行列方程式の求解を反復するため、反復回数を重ねるに従い、計算誤差が小さくなる。計算誤差0.1 kcal/molに達するには、約500分の計算量を要する。この例に見られるように、比較的少数の分子を用いて原子パラメータを決定するだけでも、必要な計算量は大きい。
【００２０】
第二の問題点は、得られた原子パラメータを溶媒和エネルギー予測に適用したときの計算精度が、事前に予測不可能なことである。一般に、計算精度は原子タイプの定義に依存すると考えられる。そのため、もし計算精度が所望の精度を下回る場合には、原子タイプを再定義し、原子パラメータ決定手順の全てをやり直す必要がある。その場合、原子パラメータ決定手順自体の反復が必要になり、計算量はさらに膨大になる。
【００２１】
この原子パラメータ決定方法は、分子を「解析対象物」、分子を構成する原子を「解析対象物を構成する要素」、原子パラメータを「要素パラメータ」、溶媒和エネルギーの実験値を単に「実験値」、溶媒和エネルギーの計算値を単に「計算値」とし、一般化可能である。
【００２２】
図６は、勾配行列を用いたパラメータ決定方法のフローチャートである。まず、入力装置601から収束閾値τ6011、解析対象物の計算に必要な情報および実験値6012、解析対象物を構成する各要素の要素タイプ6013を入力する。次に、要素パラメータ初期化ステップ602では、要素タイプごとに設定した要素パラメータを初期化する。次に、行列生成ステップ603では、全解析対象物の計算値、および全解析対象物の各要素パラメータの変化に対する計算値勾配を格納する行列を生成する。次に、行列方程式求解ステップ604では、生成した行列をA,  要素パラメータの更新値ベクトルをx、実験値ベクトルをbとしたとき、目的関数g=Σ|Ax-b|2を最小化するような更新値ベクトルxを求解する。次に、要素パラメータ更新ステップ605では、更新値ベクトルに基づき要素パラメータ値を更新する。次に、要素パラメータ計算終了判定ステップ606では、計算量の残差平方和の変化量がτ未満であるか判定する。変化量がτ以上の場合は、行列生成ステップ603、行列方程式求解ステップ604、要素パラメータ更新ステップ605を反復する。要素パラメータの変化量がτ未満の場合は、出力装置607より、収束した要素パラメータ値6071を出力する。
【００２３】
図７は、勾配行列を用いたパラメータ決定方法における、行列生成ステップ、行列方程式求解ステップ、および要素パラメータ更新ステップにおける処理の詳細である。行列生成ステップ701で実行する処理702は、全解析対象物の計算値7021、および各要素パラメータの変化に対する全解析対象物の計算値勾配7022の計算である。計算結果は行列7023に格納する。行列内で○印で示す行列要素7024の一つ一つに要する計算量が、計算値の評価一回分に相当する。行列方程式求解ステップ703で実行する処理704は、行列生成ステップ701で生成された行数が解析対象物数、列数が（要素タイプ数＋１）からなる行列7041と、（要素タイプ数＋１）個の要素数を持つ解ベクトル7042と、解析対象物数と同一の要素数を持つ実験値ベクトル7043で構成される行列方程式から定義できる目的関数gの最小化計算である。行列7041と、実験値ベクトル7043は既知であり、解ベクトル7042が求解対象の未知量である。行列方程式の求解方法は複数あるが、計算の律速段階は行列生成ステップであり、行列方程式求解ステップではないことを前提としているので、ここでは求解の具体的な方法は任意であるとする。要素パラメータ更新ステップ705で実行する処理706は、新しい要素パラメータベクトル7061の、更新前の要素パラメータベクトル7062と行列方程式求解ステップ703で得られた解ベクトルから先頭要素のみを削除したベクトル7063の和による更新である。ベクトル7063の各要素は、要素パラメータ値の次回ループ時の変化量を表している。要素パラメータを更新すると、計算値や、計算値の要素パラメータ変化に対する勾配が変化するので、再び行列生成ステップから反復する。
【００２４】
一般化されたパラメータ決定方法において、行列要素の計算が全計算の大部分を占める場合には、行列方程式の反復求解に要する計算量が膨大となる問題がある。計算量は数２で与えられる。但し、f(A(n))の関数形は解析対象物により異なる。また、パラメータ決定の結果、所望の精度を満たさない場合には、パラメータを定義し直し、パラメータ決定手順自体の反復が必要となるので、必要な計算量はさらに大きくなる。
【発明の開示】
【発明が解決しようとする課題】
【００２５】
本発明では、上記の第一の問題点を解決するため、パラメータ決定に要する計算量の削減が第一の課題となる。
【００２６】
また、本発明では、上記の第二の問題点を解決するため、所望の精度を満たすパラメータ決定に要する計算量の削減が第二の課題となる。
【課題を解決するための手段】
【００２７】
上記の第一の課題を解決するため、本発明の実施例に係る原子パラメータ決定方法では、未決定の原子パラメータ数が１以下である分子のみを選択し、選択した分子を、計算順序を制御するリストに格納し、そのリストに基づいて未決定原子パラメータを逐次決定する。そして、全ての原子パラメータが決定するまで分子の選択を反復する。
【００２８】
反復による原子パラメータ決定手順の初回では、未決定の原子パラメータ数が１である分子を選択する。そのような分子の溶媒和エネルギーの実験値を計算で再現するには、一つの原子パラメータ値を決定すればよい。したがって、原子パラメータは一意に決定できる。
【００２９】
次に、未決定の原子パラメータ数が２である分子を考える。２つの未決定パラメータのうちのどちらかが、反復の初回で決定した原子パラメータ値に基づいて仮決定できれば、未決定の原子パラメータ数は１となる。そこで、その１つの原子パラメータを決定すればよい。仮決定した原子パラメータの信頼度が高い場合には、新たに決定する原子パラメータも確度高く決定できる。
【００３０】
高い信頼度で原子パラメータを仮決定するために、各原子について、原子電荷と平均結合距離があらかじめ計算しておく。そして、原子パラメータを、原子電荷と平均結合距離を変数とする関数として表現する。原子電荷は分子と水溶媒間のクーロン相互作用強度と関連する量であり、平均結合距離は分子の結合原子数、結合次数、結合相手原子などの結合環境を反映する量である。そこで、原子電荷と平均結合距離が共に類似である原子では、原子パラメータ値も類似していると考えられる。
【００３１】
しかしながら、上記の方法で得られた原子パラメータが所望の精度、たとえば1.4 kcal/mol、で溶媒和エネルギーを決定できる保証はない。そこで、上記の第二の課題を解決するため、本発明の実施例に係る原子パラメータ決定方法では、まずシミュレーション計算で得られる溶媒和エネルギーの計算値と溶媒和エネルギーの実験値の差の、許容できる最大値である許容エネルギー誤差を設定する。そして、計算順序を制御するリストに基づいて未決定原子パラメータを逐次決定する過程において、各分子について計算値と実験値の差と定義する計算誤差の絶対値が許容誤差以内になっているかを判定する。
【００３２】
もし、計算誤差の絶対値の方が許容誤差よりも大きい場合には、計算値を実験値に近づけるため、分子に含まれる複数の原子パラメータのうち、最も溶媒和エネルギーへの寄与が大きな一つの原子パラメータに着目する。その原子パラメータの値は仮決定しているが、その値を未決定と設定し直し、新たに原子パラメータ値を決定する。この方法では、原子パラメータの数ではなく、許容誤差が固定される。
【００３３】
先に、原子パラメータの値は、原子電荷と平均結合距離を変数とする関数として表現すると述べた。上記の例に見られるように、単純な関数に基づいて仮決定した原子パラメータを、新たに未決定と設定し直し、再決定する処理を実行した場合、原子パラメータの値を表す関数はより複雑となる。
【００３４】
例として、ある原子タイプに関して二つの原子パラメータが決定しており、原子パラメータ値を原子電荷を変数とする一次関数で表現する場合を考える。同一の原子タイプに属する別の原子のパラメータ値は、その原子の原子電荷を、原子パラメータ値を表現する一次関数に代入し、仮決定できる。その仮決定した原子パラメータ値が、パラメータ決定手順の過程で未決定であると設定し直され、値が再決定された場合には、着目する原子タイプに属する３つの原子の原子パラメータ値は、原子電荷を変数とする一次関数では表現不可能となる。しかし、原子電荷と平均結合距離の両者を変数とする一次関数や、原子電荷を変数とする二種類の一次関数のようなより複雑な関数を用いれば、原子パラメータ値を表現可能である。そこで本発明では、計算誤差の絶対値が許容誤差以内となるように、原子パラメータ決定の過程で、各原子パラメータの値を表す関数を複雑化する。
【００３５】
ここまで、例として、分子の溶媒和エネルギー計算に用いる原子パラメータの決定方法を取り上げた。この例を一般化し、分子の代わりに材料、部品、構造物などの解析対象物を考えても、計算順序リストを用いて解析対象物の計算順序を制御する手段は有効であると考える。本発明の適用範囲は、計算対象が要素に分割でき、各要素に属性値が定義できる場合で、かつ、計算値や計算値のパラメータ変化に対する勾配の評価が計算の律速段階である場合である。また、計算対象について測定された特性値や物性値などの実験値とシミュレーションによって得られる計算値の差を許容誤差以内にするため、要素のパラメータ値を要素物性値の関数として表し、必要に応じて関数形を複雑化する手段も有効であると考える。この手段の適用範囲は、実験値と計算値をフィットさせるために用いるパラメータの数が増加してもよい場合である。
未決定の原子パラメータ数が１以下である分子のみを選択する手段は、高確度に決定可能な原子パラメータを優先的に決定する作用を持つ。例えば、ベンゼン分子には６つの炭素原子から構成されるが、それらは対称性の観点から等価であるので、原子パラメータ数は１となる。したがって、ベンゼンの溶媒和エネルギーの実験値を再現する原子パラメータ値は高確度に、この場合は一意に、決定可能である。
【００３６】
また、選択した分子を、計算順序を制御するリストに格納し、そのリストに基づいて一変数方程式の求解を反復し未決定原子パラメータを逐次決定する手段は、従来の技術における行列生成ステップを不要とする。そのため、この手段は、原子パラメータ決定に要する計算量を削減する作用がある。従来の技術では、全分子の溶媒和エネルギーおよび各原子パラメータの変化に対する全分子の溶媒和エネルギー勾配の計算を反復していた。一方、本発明の実施例では、分子の計算順序を決定し、各分子について溶媒和エネルギーと一つの原子パラメータの変化に対する溶媒和エネルギー勾配の計算を反復する。
【００３７】
また、計算順序を制御するリストに基づいて一変数方程式の求解を反復し未決定原子パラメータを逐次決定する過程において、ある分子について仮決定している原子パラメータ値を、未決定と設定し直し、新たに原子パラメータ値を決定する手段は、溶媒和エネルギーを所望の精度で再現する作用がある。仮決定していた原子パラメータ値を、新たに決定し直せば、溶媒和エネルギーの計算精度を向上できるからである。溶媒和エネルギーの実験値を再現するために用いるパラメータの数を増やしてもよいという条件下では、この手段を用いて原子タイプの再定義と、原子パラメータ決定手順全体のやり直しを回避できる。
【００３８】
分子の代わりに材料、部品、構造物などの解析対象物を考えた場合も、作用は同様である。未決定の要素パラメータ数が１以下である解析対象物のみを選択する手段は、高確度に決定可能な要素パラメータを優先的に決定する作用を持つ。また、選択した解析対象物を、計算順序を制御するリストに格納し、そのリストに基づいて一変数方程式の求解を反復し未決定要素パラメータを逐次決定する手段は、従来の技術における行列生成ステップを不要とする。さらに、計算順序を制御するリストに基づいて一変数方程式の求解を反復し未決定要素パラメータを逐次決定する過程において、ある解析対象物について仮決定している要素パラメータ値を、未決定と設定し直し、新たに要素パラメータ値を決定する手段は、計算値を所望の精度で再現する作用がある。
【発明の効果】
【００３９】
本発明によれば、従来法よりも計算量を小さくすることが可能である。
【発明を実施するための最良の形態】
【００４０】
以下、図面を参照して、本発明の実施の形態について説明する。
（実施例１）
図８、および図９を用いて、本発明の原子パラメータ決定方法を説明する。
【００４１】
図８は、原子パラメータ決定方法のフローチャートである。まず、入力装置801から収束エネルギー閾値τ8011、分子の計算または実験により決定した三次元座標情報と溶媒和エネルギーの実験値8012、分子構成する各原子の原子タイプ、および、各原子の原子電荷と平均結合距離8013を入力する。次に、初期化ステップ802では、原子タイプごとに設定した原子パラメータの初期化ステップ8021と、全分子を未決定原子パラメータ数の昇順に並べ替えるステップ8022を実行する。計算順序リスト生成ステップ803では、未決定の原子パラメータを含む分子から、未決定の原子パラメータ数が１である分子を選択するステップ8031と、選択された分子を原子パラメータ数の昇順に並べ替えた結果を格納する計算順序リストを生成するステップ8032を実行する。次に、一変数方程式求解ステップ804では、計算順序リストの最上位分子について、溶媒和エネルギー計算値の変化量がτ未満となるまで、一変数方程式の求解を反復し、未決定原子パラメータを決定するステップ8041と、計算順序リストの更新ステップ8042を実行する。計算順序リスト終了判定ステップ8043で、計算順序リストに最上位の分子が存在する場合にはステップ8041に戻る。存在しない場合には、原子パラメータ更新ステップ805へと進む。
【００４２】
次に、原子パラメータ更新ステップ805では、決定した原子パラメータに対応する原子タイプに関して、原子電荷と平均結合距離の関数である原子パラメータの表現式を更新するステップ8051と、未決定分子を未決定原子パラメータの昇順に並べ替えるステップ8052と、未決定原子パラメータ数が同一の分子を原子パラメータ数の降順に並べ替えるステップ8053と、未決定分子を更新するステップ8054を実行する。原子パラメータ計算終了判定ステップ806では、未決定分子が存在するか判定する。存在する場合には、計算順序リスト生成ステップ803、一変数方程式求解ステップ804、原子パラメータ更新ステップ805を繰り返す。存在しない場合には、計算誤差算出ステップ807で、得られた原子パラメータを用いて各分子の計算誤差を算出し、出力装置808より、決定した原子パラメータ値8081を出力する。
【００４３】
図９は、原子パラメータ決定方法における、計算順序リスト生成ステップ、一変数方程式求解ステップ、および原子パラメータ更新ステップにおける処理の詳細である。例として、分子数が９、原子タイプ数が８の場合における各ステップの一回目と二回目の繰り返し過程を説明する。計算順序リスト生成ステップ901で実行する一回目の処理902では、あらかじめ未決定原子パラメータ数9021の昇順に並べられている分子から、未決定原子パラメータ数が１である分子9022（分子８、分子３、分子４）を選択する。これらの分子の計算順序は、分子の並び順と同じ分子８、分子３、分子４の順であり、この順序を計算順序リスト9023に格納する。一変数方程式求解ステップ903で実行する一回目の処理904では、計算順序リストに基づき未決定原子パラメータを決定する。最初の計算は分子８であるので、まず、原子パラメータra,8を用いて溶媒和エネルギー計算をし、実験値E8を再現するようにra,8の値を調整する。この過程では、ra,8を変化させつつ実験値を再現できる値を探索するので、一変数方程式を反復求解が必要となる。選択された残り二つの分子についても同様の処理を実行する。原子パラメータ更新ステップ905で実行する一回目の処理906では、決定した原子パラメータが対応する原子タイプ9061であるa, b, およびhに関して、原子電荷と平均結合距離の関数である原子パラメータの表現式を更新する。原子パラメータの表現式は、原子電荷と平均結合距離のどちらか一方の関数であってもよい。例えば、原子タイプaの原子パラメータ9062はra(q)=c0+c1qと表現できる。ここで、c0,c1は係数である。次に、原子タイプaに分類されている未決定の原子パラメータ値を、ra(q)=c0+c1qに原子電荷qを代入して得る。ただし、原子電荷がパラメータ決定に用いられた原子の原子電荷とある閾値以上離れた原子9063に関しては、原子パラメータ値を設定しない選択肢もある。値が割り当てられている原子910を●、値が割り当てられていない、すなわち未決定原子911を○で示す。次に、この時点で未決定の分子である、分子５、分子６、分子７、分子２、分子１、分子９の６つを未決定原子パラメータ数9064の昇順に並べ替える。また、未決定原子パラメータ数が等しい場合には、原子パラメータ数9065の降順に並べ替える。このとき、未決定パラメータ数が１以上の分子が、新たな未決定分子となる。分子２は、この時点で未決定原子パラメータ数が０となるので、原子パラメータの決定には用いない。計算順序リスト生成ステップ901で実行する二回目の処理907では、あらかじめ並べられている分子から、未決定原子パラメータ数9071が１である分子9072（分子５、分子６）を選択する。これらの分子の計算順序は分子５、分子６の順であり、この順序を計算順序リスト9073に格納する。一変数方程式求解ステップ903で実行する二回目の処理908では、計算順序リストに基づいた未決定原子パラメータを決定する。まず、値が固定された原子パラメータra,5と、未決定原子パラメータrd,5を用いて溶媒和エネルギー計算をし、実験値E5を再現するようにrd,5を決定する。この過程では、rd,5のみを変化させつつ実験値を再現できる値を探索するので、二つの原子パラメータから構成される分子であっても、実行する処理は一変数方程式を反復求解である。選択された残り一つ��分子についても同様の処理を実行する。原子パラメータ更新ステップ905で実行する二回目の処理909では、決定した原子パラメータが対応する原子タイプ9091であるa, b, h, dおよびeに関する、原子電荷と平均結合距離の関数である原子パラメータの表現式を更新する。例えば、原子タイプaの原子パラメータ5092はra(q)=c0+c1qから、ra(q)=c'0+c'1qへと更新する。表現式が更新した後に、該当する原子タイプに分類される原子の原子パラメータ値を決定し、未決定の分子を並べ替える操作を反復の一回目同様に実行する。この処理を未決定の分子がなくなるまで反復する。
【００４４】
以下に、図８に示した原子パラメータ決定方法のフローチャートの一部のステップについて、処理の実施例を具体的に説明する。
【００４５】
分子の三次元座標情報および溶媒和エネルギーの実験値の入力8012の例を図１０に示す。また、分子を構成する各原子の原子タイプtおよび原子電荷q、平均結合距離lの入力8013の例を図１１に示す。これらの入力は相互に関連しているため、図１０と図１１をまとめて説明する。まず、入力情報は、分子に帰属する情報と、原子に帰属する情報に分類できる。分子に帰属する入力情報として分子の三次元座標情報1001がある。例としてパラキシレン分子の構造1002を示した。実際の入力は分子を構成する各原子の三次元座標である。また、分子に帰属する入力情報として溶媒和エネルギー実験値1003がある。パラキシレン分子の溶媒和エネルギーの実験値1004は-0.80 kcal/molである。決定すべき原子パラメータ数1005は、パラキシレン分子の例1006を用いると、白丸、黒丸、斜線丸の３種類となる。同一の印が付与された原子は対称性により等価であるので、同一の原子パラメータ値を持つ。この情報はパラメータ決定計算に必須ではないが、あらかじめ用意しておくと便利である。原子に帰属する情報としては、原子タイプ1101、原子電荷1102、平均結合距離1103がある。パラキシレン分子には、三種類の原子パラメータがあり、それぞれについて原子タイプ、原子電荷、平均結合距離を入力する。具体的には、原子電荷は、分子の静電ポテンシャルを最も良く近似するように、分子を構成する各原子核上に電荷を割り付けるアルゴリズムによって決定した。平均結合距離は、各原子が直接結合している原子に対する結合距離の平均値として計算した。
【００４６】
原子タイプごとに設定した原子パラメータrt(q,l)の初期化8021の例を図１２に示す。ここでも、取り扱う情報を分子に帰属する情報と、原子に帰属する情報に分類すると便利である。原子に帰属する情報1201には、各原子パラメータの値1202と、原子パラメータの状態1203がある。原子パラメータの状態は、値あり1204、値なし1205のいずれかである。初期化値としては、全ての原子パラメータに関して値は定義されておらず、状態は「値なし」とする。分子に帰属する情報1206は、原子に帰属する情報の集計により算出できる。分子に含まれる「値あり」の原子パラメータ数1207は原子パラメータの状態が「値あり」である原子パラメータ数であり、「値なし」の原子パラメータ数1208も同様に計算できる。分子の計算状態1209は「未」か「完」のいずれかであり、全ての原子パラメータが「値あり」の場合、その分子の計算状態1209は「完」となる。分子に帰属する情報の初期値1210は、原子に帰属する情報の初期値の集計結果に過ぎない。「値あり」の原子パラメータ数の初期値1211は０であり、「値なし」の原子パラメータの数1212はパラキシレンの場合、原子パラメータの数に等しい３である。計算状態の初期値1213は「未」である。
【００４７】
本発明の実施例は、計算順序リストを用いて分子の溶媒和エネルギー計算の順序を制御する部分に特徴がある。そこで、一変数方程式求解ステップ804におけるハードウェアの制御方法を、図１３を用いて説明する。
【００４８】
計算順序リスト1301の最上位分子について以下の処理を行う。まず、入力ファイル生成ステップ1302で、あらかじめ入力されている分子の三次元座標情報と、「値あり」の原子パラメータ情報に基づいて、溶媒和エネルギー計算用の入力ファイル1303を作成する。「値なし」の原子パラメータには、適当な値を設定する必要がある。例えば、元素ごとに定められているヴァンデルワールス半径値を用いてもよい。溶媒和エネルギー計算ジョブ投入ステップ1304では、計算機1305にインストールされた溶媒和エネルギー計算プログラム1306に、入力ファイル1303を渡し、溶媒和エネルギー計算を実行する。計算結果出力ステップ1307では、溶媒和エネルギー計算プログラム1306により計算結果が出力ファイル1308に書き出す。出力結果処理ステップ1309では、あらかじめ入力されている分子の溶媒和エネルギー実験値と計算値の比較に基づき、原子パラメータを更新する。入力ファイル生成ステップ1302から出力結果処理ステップ1309までを、溶媒和エネルギー計算値の変化量があらかじめ入力された収束エネルギー閾値未満となるまで反復する。反復終了の時点で、計算順序リスト最上位分子を削除し、計算順序第二位以下の分子を順次繰り上げる。上記の一連の処理を計算順序リストが空になるまで続行する。
【００４９】
決定した原子パラメータに対応する原子タイプの原子パラメータrt(q,l)更新ステップ8051の実施例を図１４、図１５を用いて説明する。
【００５０】
図１４には、原子パラメータ値1401を、原子電荷のような、一つの原子属性値1402を変数とする一次関数で表現する場合の原子パラメータの更新例を示す。図中黒丸は原子パラメータが「値あり」である原子1403を示し、白丸は原子パラメータが「値なし」である原子1404を示す。まず、「値あり」の原子1403を用いて原子パラメータ値と原子属性値を関連付ける一次関数1405を定義する。次に、「値なし」の原子の原子属性値1406を一次関数に代入し、原子パラメータ値を割り当てる。この時、「値あり」である原子の原子属性値と大きくかけ離れた原子属性値を持つ「値なし」の原子1407に割り当てられた値は必ずしも高精度ではない可能性がある。そのような原子の原子パラメータは、その原子が属する分子の最後の原子パラメータとして決定する。
【００５１】
図１５には、原子パラメータ値1501を、原子電荷1502と平均結合距離1503のような、二つの原子属性値を変数とする一次関数で表現する場合の原子パラメータの更新例を示す。図中黒丸(および灰色丸)は原子パラメータが「値あり」である原子1504を示し、白丸は原子パラメータが「値なし」である原子1505を示す。まず、「値あり」の原子1504を用いて原子パラメータ値と原子属性値を関連付ける一次関数1506を定義する。次に、「値なし」の原子の原子電荷と平均結合距離1507を一次関数に代入し、原子パラメータ値を割り当てる。この時、「値あり」である原子の原子属性値と大きくかけ離れた原子属性値を持つ「値なし」の原子1508に割り当てられた値は必ずしも高精度ではない可能性がある。そのような原子の原子パラメータは、その原子が属する分子の最後の原子パラメータとして決定する。
【００５２】
最後に従来の原子パラメータの決定方法と、本発明による原子パラメータ決定方法における計算量の比較を行う。従来法における計算量は式１に記載している。本発明では、分子を順次選択し、各分子について溶媒和エネルギー計算値の変化量が収束エネルギー閾値未満となるまで計算を反復する。分子ごとに異なる反復回数をI(n)とする。但し、原子パラメータ決定に用いない分子も存在する。そこで、原子パラメータ決定に用いる分子では１、用いない分子では０と定義するδ関数を用いると、計算量は（数７）で計算できる。
【００５３】
【数７】
2007323099.tif 000009
【００５４】
数２と数７を比較すると、三つの点で、本発明の実施例の方法は従来法よりも計算量が少ないことがわかる。第一に、従来法では全ての分子の計算をするのに対し、本発明の実施例ではδ関数の存在により、一部の分子の計算のみを行うため計算量が少ない。第二に、従来法では数２の２行目右辺における括弧内の第一項が”I”であるのに対し、本発明の実施例では数７右辺における括弧内の第一項が”1”であるため計算量が少ない。第三に、従来法では括弧内の第二項が”I・p(n)”であるのに対し、本発明の実施例では括弧内の第一項が”I(n)”である。多変量非線形問題である行列方程式の求解に要する反復回数Iは、一変量非線形問題である一変数方程式の求解に要する反復回数I(n)よりも大きいと考えられるため、本発明の実施例の方が計算量が少ないと言える。
（実施例２）
図１６〜図１８を用いて１０分子を用いた、６つの原子パラメータの決定例を説明する。
【００５５】
図１６に、本実施例の適用に用いた原子タイプ1601、原子パラメータ値の決定過程1602、および各分子が含有する原子タイプ1603を示す。原子タイプ1601としてCH3, CH2, NH2の３種類を定義した。但し、CH3,とCH2の原子パラメータは平均結合距離lを変数とする一次関数で表現し、NH2の原子パラメータは原子電荷qを変数とする一次関数で表現する。各一次関数は変数を２つ含むので、決定すべき変数の数は６となる。原子パラメータ値は初期段階では未決定である。以下に原子パラメータの決定過程を説明する。
【００５６】
まず図１６Aにおいて、未決定原子パラメータ数が１である分子1604（番号１，９，４，１０，３の分子）を選択した。そして選択した各分子について、溶媒和エネルギーの実験値を再現する原子パラメータを決定した。
【００５７】
次に図１６Bにおいて、未決定の原子パラメータに対して値を割り当てた。原子タイプCH3に関しては、番号１と番号９の分子について得た原子パラメータに基づき、原子パラメータrCH31605をrCH3=-2.0323*l+4.5674と決定した。同様に、原子パラメータrCH31606をrCH2=3.8382*l-2.9689と決定した。原子タイプNH2に関しては、この時点では、原子パラメータを一次関数で表現するのに必要な二つの変数が決まらないため、「未決定」のままとした。原子タイプCH3と、CH2に関しては、決定した一次関数に平均結合距離lを代入し、未決定の原子パラメータ1607に値を割り当てた。その結果、未決定原子パラメータ数1608が分子７、分子５に関しては０となり、分子２、分子６、分子８に関しては１となった。
【００５８】
次に、図１６Cにおいて、未決定原子パラメータ数が１である分子1609（番号２，６，８の分子）を選択した。そして選択した各分子について、溶媒和エネルギーの実験値を再現する原子パラメータを決定した。
【００５９】
最後に、図１６Dにおいて、原子タイプNH2に関して決定した４つの原子パラメータ1610と平均結合距離から、原子パラメータrNH21611をrNH2=0.2943*q+1.8561と決定した。番号３、２、６、８の各分子については、各分子について決定した原子パラメータと、一次関数に原子電荷を代入して得られる原子パラメータが異なるので、最後に後者により得られた原子パラメータ値を用いて、計算誤差を算出した。
【００６０】
図１７に、上記「従来の技術」の欄で説明した従来法と本実施例の計算量および計算誤差の比較を示す。図中横軸は一台のコンピュータを用いた原子パラメータ決定に要した計算時間1701と定義する計算量、縦軸は各分子の絶対計算誤差の平均値と定義する計算誤差1702である。ひし形プロット1703は従来法における計算量と計算誤差の関係であり、長方形プロット1704は本実施例における計算量と計算誤差の関係である。本実施例では、分子を順次選択し、全ての原子パラメータが決定した時点で計算を終了する。要した計算量は57分であり、約0.1 kcal/molの計算誤差を得た。したがって、本実施例の方法は、本実施例における計算条件下では、計算誤差0.1 kcal/molを達成するのに要する計算量が約９倍少ないと言える。
【００６１】
図１８に溶媒和エネルギー予測精度の比較を示す。６分子に対して、従来法による１２回反復後の原子パラメータと、本実施例の方法で得られた原子パラメータをそれぞれ適用した。選択した分子1801は、原子パラメータ決定に用いた１０分子と類似構造を持つ。図には、実験値1802、従来法による計算誤差1803、本実施例による計算誤差1804、および両方法における平均絶対誤差1805を示す。本実施例の方法は、原子パラメータを常に１つずつ決定するので、従来法による原子パラメータ決定方法に対する近似と考えてもよい。だが、得られた平均絶対誤差1805は、本実施例の方が従来法よりも僅かに小さい。したがって、原子パラメータ決定に要する計算量を削減する本実施例は、溶媒和エネルギーの予測精度も保持すると言える。むしろ、原子パラメータを原子電荷や平均結合距離の関数として表現する手法は、高精度な溶媒和エネルギー計算に有効であると考えられる。
（実施例３）
本発明の本実施例に係る原子パラメータ決定方法は、分子を「解析対象物」、分子を構成する原子を「解析対象物を構成する要素」、原子パラメータを「要素パラメータ」、溶媒和エネルギーの実験値を単に「実験値」、溶媒和エネルギーの計算値を単に「計算値」とし、一般化可能である。以下に、図１９、図２０を用いて、一般化されたパラメータ決定方法を、従来の方法と対比し、説明する。
【００６２】
図１９は、本実施例のパラメータ決定方法のフローチャートである。まず、入力装置1901から収束閾値τ19011、解析対象物の計算に必要な情報および実験値19012、解析対象物構成する各要素の要素タイプ、および、各要素の要素属性値19013を入力する。次に、初期化ステップ1902では、要素タイプごとに設定した要素パラメータの初期化ステップ19021と、全解析対象物を未決定要素パラメータ数の昇順に並べ替えるステップ19022を実行する。計算順序リスト生成ステップ1903では、未決定の要素パラメータを含む解析対象物から、未決定の要素パラメータ数が１である解析対象物を選択するステップ19031と、選択解析対象物の並び順を計算順序として決定し、計算順序リストを生成するステップ19032を実行する。次に、一変数方程式求解ステップ1904では、計算順序リストの最上位解析対象物について、計算値の変化量がτ未満となるまで、一変数方程式の求解を反復し、未決定要素パラメータを決定するステップ19041と、計算順序リストの更新ステップ19042を実行する。計算順序リスト終了判定ステップ19043で、計算順序リストに最上位の解析対象物が存在する場合にはステップ19041に戻る。存在しない場合には、要素パラメータ更新ステップ1905へと進む。次に、要素パラメータ更新ステップ1905では、決定した要素パラメータに対応する要素タイプに関して、要素属性値の関数である要素パラメータの表現式を更新するステップ19051と、未決定解析対象物を未決定要素パラメータの昇順に並べ替えるステップ19052と、未決定要素パラメータ数が同一の解析対象物を要素パラメータ数の降順に並べ替えるステップ19053と、未決定解析対象物を更新するステップ19054を実行する。要素パラメータ計算終了判定ステップ1906では、未決定解析対象物が存在するか判定する。存在する場合には、計算順序リスト生成ステップ1903、一変数方程式求解ステップ1904、要素パラメータ更新ステップ1905を繰り返す。存在しない場合には、計算誤差算出ステップ1907で、得られた要素パラメータを用いて各解析対象物の計算誤差を算出し、出力装置1908より、決定した要素パラメータ値19081を出力する。
【００６３】
図２０は、本実施例のパラメータ決定方法における、計算順序リスト生成ステップ、一変数方程式求解ステップ、および要素パラメータ更新ステップにおける処理の詳細である。例として、解析対象物数が９、要素タイプ数が８の場合における各ステップの一回目と二回目の繰り返し過程を説明する。計算順序リスト生成ステップ2001で実行する一回目の処理2002では、あらかじめ未決定要素パラメータ数20021の昇順に並べられている解析対象物から、未決定要素パラメータ数が１である解析対象物20022（解析対象物８、解析対象物３、解析対象物４）を選択する。これらの解析対象物の計算順序は、解析対象物の並び順と同じ解析対象物８、解析対象物３、解析対象物４の順であり、この順序を計算順序リスト20023に格納する。一変数方程式求解ステップ2003で実行する一回目の処理2004では、計算順序リストに基づき未決定要素パラメータを決定する。最初の計算は解析対象物８であるので、まず、要素パラメータra,8を用いて計算値を取得し、実験値E8を再現できるようにra,8を決定する。この過程では、ra,8を変化させつつ実験値を再現できる値を探索するので、一変数方程式を反復求解が必要となる。選択された残り二つの解析対象物についても同様の処理を実行する。要素パラメータ更新ステップ2005で実行する一回目の処理2006では、決定した要素パラメータが対応する要素タイプ20061であるa, b, およびhに関して、要素属性値の関数である要素パラメータの表現式を更新する。例えば、要素タイプaの要素パラメータ20062はra(q)=c0+c1p1と表現できる。ここで、c0,c1は係数である。次に、要素タイプaに分類されている未決定の要素パラメータ値を、ra(q)=c0+c1p1に要素属性値p1を代入して得る。ただし、要素属性値がパラメータ決定に用いられた要素の要素属性値とある閾値以上離れた要素20063に関しては、要素パラメータ値を設定しない選択肢もある。値が割り当てられている要素2010を●、値が割り当てられていない、すなわち未決定要素2011を○で示す。次に、この時点で未決定の解析対象物である、解析対象物５、解析対象物６、解析対象物７、解析対象物２、解析対象物１、解析対象物９の６つを未決定要素パラメータ数20064の昇順に並べ替える。また、未決定要素パラメータ数が等しい場合には、要素パラメータ数20065の降順に並べ替える。このとき、未決定パラメータ数が１以上の解析対象物が、新たな未決定解析対象物となる。解析対象物２は、この時点で未決定要素パラメータ数が０となるので、要素パラメータの決定には用いない。計算順序リスト生成ステップ2001で実行する二回目の処理2007では、あらかじめ並べられている解析対象物から、未決定要素パラメータ数20071が１である解析対象物20072（解析対象物５、解析対象物６）を選択する。これらの解析対象物の計算順序は解析対象物５、解析対象物６の順であり、この順序を計算順序リスト20073に格納する。一変数方程式求解ステップ2003で実行する二回目の処理2008では、計算順序リストに基づいた未決定要素パラメータを決定する。まず、値が固定された要素パラメータra,5と、未決定要素パラメータrd,5を用いて計算値を取得し、実験値E5を再現するようにrd,5を決定する。この過程では、rd,5のみを変化させつつ実験値を再現できる値を探索するので、二つの要素パラメータから構成される解析対象物であっても、実行する処理は一変数方程式を反復求解である。選択された残り一つの解析対象物についても同様の処理を実行する。要素パラメータ更新ステップ2005で実行する二回目の処理2009では、決定した要素パラメータが対応する要素タイプ20091であるa, b, h, dおよびeに関する、要素属性値の関数である要素パラメータの表現式を更新する。例えば、要素タイプaの要素パラメータ20092はra(q)=c0+c1p1から、ra(q)=c'0+c'1p1へと更新する。表現式が更新した後に、該当する要素タイプに分類される要素の要素パラメータ値を決定し、未決定の解析対象物を並べ替える操作を反復の一回目同様に実行する。この処理を未決定の解析対象物がなくなるまで反復する。
【００６４】
最後に従来のパラメータの決定方法と、本実施例によるパラメータ決定方法における計算量の比較を行う。従来法における計算量は式１に記載している。本実施例では、解析対象物を順次選択し、各解析対象物について計算値の変化量が収束閾値未満となるまで計算を反復する。解析対象物ごとに異なる反復回数をI(n)とする。但し、要素パラメータ決定に用いない解析対象物も存在する。そこで、要素パラメータ決定に用いる解析対象物では１、用いない解析対象物では０と定義するδ関数を用いると、計算量は（数８）、
【００６５】
【数８】
2007323099.tif 000010
となる。
【００６６】
数２と数８を比較すると、三つの点で、本実施例の方法は従来法よりも計算量が少ないことがわかる。第一に、従来法では全ての解析対象物の計算をするのに対し、本実施例ではδ関数の存在により、一部の解析対象物の計算のみを行うため計算量が少ない。第二に、従来法では数２の２行目右辺における括弧内の第一項が”I”であるのに対し、本実施例では数８右辺における括弧内の第一項が”1”であるため計算量が少ない。第三に、従来法では括弧内の第二項が”I・p(n)”であるのに対し、本実施例では括弧内の第一項が”I(n)”である。多変量非線形問題である行列方程式の求解に要する反復回数Iは、一変量非線形問題である一変数方程式の求解に要する反復回数I(n)よりも大きいと考えられるため、本実施例の方が計算量が少ないと言える。
（実施例４）
図２１、および図２２を用いて、本発明の本実施例に係る許容誤差固定による原子パラメータ決定方法を説明する。
【００６７】
図２１は、許容誤差固定による原子パラメータ決定方法のフローチャートである。まず、入力装置2101から許容エネルギー誤差ε21011、収束エネルギー閾値τ21012、分子の三次元座標情報および溶媒和エネルギーの実験値21013、分子構成する各原子の原子タイプ、および、各原子の原子電荷と平均結合距離21014を入力する。次に、初期化ステップ2102では、原子タイプごとに設定した原子パラメータの初期化ステップ21021と、全分子を未決定原子パラメータ数の昇順に並べ替えるステップ21022を実行する。計算順序リスト生成ステップ2103では、未決定の原子パラメータを含む分子から、未決定の原子パラメータ数が１以下である分子を選択するステップ21031と、選択分子の並び順を計算順序として決定し、計算順序リストを生成するステップ21032を実行する。次に、一変数方程式求解ステップ2104では、計算順序リストの最上位分子が未決定原子パラメータを含むかの判定ステップ21041で含むと判定された場合、溶媒和エネルギー計算値の変化量がτ未満となるまで、一変数方程式の求解を反復し、未決定原子パラメータを決定するステップ21042を実行する。つぎに、判定ステップ21041の結果に関わらず、決定した原子パラメータを用いて該分子の計算値と実験値の差である計算誤差を算出するステップ21043と、計算誤差の絶対値と許容エネルギー誤差を比較するステップ21404と、該計算誤差の方が許容エネルギー誤差を超える場合には、仮決定の原子パラメータの１つを未決定の原子パラメータとして新規に設定するステップ21045を経てステップ21041に戻り、計算誤差の方が許容エネルギー誤差以下の場合には計算順序リスト更新ステップ21046を実行する。計算順序リスト終了判定ステップ21047で、計算順序リストに最上位の分子が存在する場合にはステップ21041に戻る。存在しない場合には、原子パラメータ更新ステップ2105へと進む。次に、原子パラメータ更新ステップ2105では、決定した原子パラメータに対応する原子タイプに関して、原子電荷と平均結合距離の関数である原子パラメータの表現式を更新するステップ21051と、未決定の原子の一部または全部の原子パラメータ値を仮決定するステップ21052と、未決定分子を未決定原子パラメータの昇順に並べ替えるステップ21053と、未決定原子パラメータ数が同一の分子を原子パラメータ数の降順に並べ替えるステップ21054と、未決定分子を更新するステップ21055を実行する。原子パラメータ計算終了判定ステップ2106では、未決定分子が存在するか判定する。存在する場合には、計算順序リスト生成ステップ2103、一変数方程式求解ステップ2104、原子パラメータ更新ステップ2105を繰り返す。存在しない場合には、出力装置2107より、決定した原子パラメータ値21071を出力する。
【００６８】
図２２は、原子パラメータ決定方法における、計算順序リスト生成ステップ、一変数方程式求解ステップ、および原子パラメータ更新ステップにおける処理の詳細である。例として、分子数が９、原子タイプ数が８の場合における各ステップの一回目と二回目の繰り返し過程を説明する。計算順序リスト生成ステップ2201で実行する一回目の処理2202では、あらかじめ未決定原子パラメータ数22021の昇順に並べられている分子から、未決定原子パラメータ数が１以下である分子22022（分子８、分子３、分子４）を選択する。これらの分子の計算順序は、分子の並び順と同じ分子８、分子３、分子４の順であり、この順序を計算順序リスト22023に格納する。一変数方程式求解ステップ2203で実行する一回目の処理2204では、計算順序リストに基づき未決定原子パラメータを決定する。最初の計算は分子８であるので、原子パラメータra,8を用いて溶媒和エネルギー計算をし、実験値E8を再現できるようにra,8を決定する。この過程では、ra,8を変化させつつ実験値を再現できる値を探索するので、一変数方程式の反復求解が必要となる。選択された残り二つの分子についても同様の処理を実行する。原子パラメータ更新ステップ2205で実行する一回目の処理2106では、決定した原子パラメータが対応する原子タイプ22061であるa, b, およびhに関して、原子電荷と平均結合距離の関数である原子パラメータの表現式を更新する。原子パラメータの表現式は、原子属性値を変数に持たない関数であってもよい。例えば、原子タイプaの原子パラメータ22062はra(q,l)=c0=ra,8と表現できる。原子タイプaに分類されている原子の原子パラメータ値は、ra(q,l)=c0と仮決定できる。ただし、原子電荷あるいは平均結合距離がパラメータ決定に用いられた原子の原子電荷や平均結合距離とある閾値以上離れた原子22063に関しては、原子パラメータ値を仮決定しない選択肢もある。原子パラメータの値が決定している原子2210を●、仮決定している原子2211を▲、値が割り当てられていない、すなわち未決定原子2212を○で示す。次に、この時点で未決定の分子である、分子５、分子６、分子７、分子２、分子１、分子９の６つを未決定原子パラメータ数22064の昇順に並べ替える。また、未決定原子パラメータ数が等しい場合には、原子パラメータ数22065の降順に並べ替える。これら６分子全てが、新たな未決定分子となる。計算順序リスト生成ステップ2201で実行する二回目の処理2207では、あらかじめ並べられている分子から、未決定原子パラメータ数22071が１以下である分子22072（分子２、分子５、分子６）を選択する。これらの分子の計算順序は分子２、分子５、分子６の順であり、この順序を計算順序リスト22073に格納する。一変数方程式求解ステップ2203で実行する二回目の処理2208では、計算順序リストに基づいた未決定原子パラメータを決定する。未決定パラメータが存在しない分子２の場合について特に説明する。まず、与えられた原子パラメータを用いて溶媒和エネルギーを計算する。次に、計算値と実験値の差と定義する計算誤差の絶対値が許容エネルギー誤差未満となっているかを確認する。計算誤差の絶対値が許容エネルギー誤差以上の場合は、分子２に含まれる二つ原子パラメータra,2とrh,2のうちいずれかを未決定原子パラメータとして新規に設定する。選択する��子パラメータは、両原子パラメータの変化に対する溶媒和エネルギーの勾配が大きい方とする。選択されたパラメータがra,2である場合には、溶媒和エネルギーの実験値E2を再現するようにra,2の値を決定する。原子パラメータ更新ステップ2205で実行する二回目の処理2209では、決定した原子パラメータが対応する原子タイプ22091であるa, b, h, dおよびeに関する、原子電荷と平均結合距離の関数である原子パラメータの表現式を更新する。例えば、原子タイプaの原子パラメータ22092はra(q,l)=c0=ra,8であったが、分子２の計算結果から、原子タイプaの原子パラメータには新たに決定したra,2の値も用いる必要があるので、より複雑な関数形を用いる必要がある。例えば、ra(q,l)=c'0+c'1qと更新してもよい。表現式を更新した後に、該当する原子タイプに分類される原子の原子パラメータ値を決定し、未決定の分子を並べ替える操作を反復の一回目同様に実行する。この処理を未決定の分子がなくなるまで反復する。
【００６９】
以下に、図２１に示した許容誤差固定による原子パラメータ決定方法のフローチャートのステップのうち、実施例１における処理と異なる部分について、処理の実装例を具体的に説明する。
【００７０】
図２３に、計算順序リスト更新ステップの詳細を示す。開始2301時点では、計算順序リスト最上位の分子の計算が終了している。未決定パラメータの新規設定判定ステップ2302では、該分子の未決定パラメータ決定過程において、仮決定の原子パラメータの一つを未決定原子パラメータとして新規設定する処理が行われたかを判定する。判断結果がyesの場合は、計算順序リストからすべての分子を削除する処理2303を実行する。判断結果がnoの場合は、計算順序リストから最上位分子を削除し、残った分子を順次繰り上げる処理2304を実行する。次に計算順序リスト終了判定ステップ2305では、計算順序リストが空であるかを判定する。空である場合2306には、原子パラメータ更新ステップへ進む。空でない場合2307には、未決定原子パラメータ決定ステップへ戻る。
【００７１】
決定した原子パラメータに対応する原子タイプの原子パラメータrt(q,l)更新ステップ21051の実施例を図２４、図２５を用いて説明する。
【００７２】
図２４Aには、原子パラメータ値2401を、原子電荷のような、一つの原子属性値2402を変数とする関数で表現する場合の原子パラメータの更新例を示す。図中黒丸は原子パラメータの「値あり」の原子2403を示し、白丸は原子パラメータの「値なし」の原子2404を示す。まず、「値あり」の原子2403を用いて原子パラメータ値と原子属性値を関連付ける一次関数2405を定義する。次に、「値なし」の原子の原子属性値2406を一次関数に代入し、原子パラメータ値を仮決定する。
【００７３】
図２４Bには、仮決定した原子パラメータ値を用いて、計算順序リスト生成ステップ、一変数方程式求解ステップ、原子パラメータ更新ステップを実行したあとの結果の例を示す。仮決定した原子パラメータ値がそのまま決定値となった原子2407と、仮決定した原子パラメータ値が未決定原子パラメータとして設定され、原子パラメータが再決定した原子2408が存在する。これらの原子の原子パラメータ値を原子属性値の一次関数で記述すると、すでに決定した原子パラメータの値が変化し、原子パラメータ決定済みの分子の溶媒和エネルギーもが変化してしまう。しかし、原子属性値に基づいて三つの区域2409を定義し、それぞれの区域に対して一次関数2410を定義すれば、原子パラメータ決定済みの分子の溶媒和エネルギーを変えずに、新規に決定した原子パラメータを関数として表現できる。
【００７４】
図２５Aには、原子パラメータ値2501を、原子電荷2502と平均結合距離2503のような、二つの原子属性値を変数とする関数で表現する場合の原子パラメータの更新例を示す。図中黒丸は原子パラメータの「値あり」の原子2504を示し、白丸は原子パラメータの「値なし」の原子2505を示す。まず、「値あり」の原子2504を用いて原子パラメータ値と原子属性値を関連付ける一次関数2506を定義する。次に、「値なし」の原子の原子電荷と平均結合距離2507を一次関数に代入し、原子パラメータ値を仮決定する。
【００７５】
図２５Bには、仮決定した原子パラメータ値を用いて、計算順序リスト生成ステップ、一変数方程式求解ステップ、原子パラメータ更新ステップを実行したあとの結果の例を示す。仮決定した値がそのまま決定値となった原子2508と、仮決定した値が未決定原子パラメータとして設定され、原子パラメータが再決定した原子2509が存在する。これらの原子の原子パラメータ値を原子属性値の一次関数で記述すると、すでに決定した原子パラメータの値が変化し、原子パラメータ決定済みの分子の溶媒和エネルギーも変化してしまう。しかし、原子属性値に基づいて、原子電荷と平均結合距離が成す平面上に二つの区域2510を定義し、それぞれの区域に対して一次関数2511を定義すれば、原子パラメータ決定済みの分子の溶媒和エネルギーを変えずに、新規に決定した原子パラメータを関数として表現できる。
【００７６】
最後に従来の原子パラメータの決定方法と、本実施例による原子パラメータ決定方法における計算量の比較を行う。従来法における計算量は数２に記載している。本実施例では、分子を順次選択し、各分子について溶媒和エネルギー計算値の変化量が収束エネルギー閾値未満となるまで計算を反復する。分子ごとに異なる反復回数をI(n)とする。但し、未決定原子パラメータ数が０である分子も存在する。そこで、未決定原子パラメータ決定処理を必要とする分子では１、必要としない分子では０と定義するδ関数を用いると、計算量は（数９）
【００７７】
【数９】
2007323099.tif 000011
となる。
【００７８】
数２と数９を比較すると、二つの点で、本実施例の方法は従来法よりも計算量が少ないことがわかる。第一に、従来法では式１の２行目右辺における括弧内の第一項が”I”であるのに対し、本実施例では括弧内の第一項が”1”であるため計算量が少ない。第二に、従来法では括弧内の第二項が”I・p(n)”であるのに対し、本実施例では括弧内の第一項が”I(n)・δ(n)”であるため以下の二つの理由により、計算量が少ない。一つ目に、多変量非線形問題である行列方程式の求解に要する反復回数Iは、一変量非線形問題である一変数方程式の求解に要する反復回数I(n)よりも大きいと考えられる。二つ目に、分子nに含まれる原子タイプ数p(n)の最低値が１であるのに対し、δ関数は０か１の値である。
【００７９】
さらに、従来法と本実施例の方法の計算量の差は、分子数Nが大きい場合に、増大する。以下にその理由を説明する。
【００８０】
本実施例の方法では、原子パラメータを原子電荷と平均結合距離を変数とする関数で表現する。その関数形は、分子を選択し、未決定のパラメータを決定する過程で複雑化する。しかし、原子がとりうる原子電荷や平均結合距離の範囲は有限である。そのため、関数形はある一定段階まで複雑化した後は、あらゆる原子電荷と平均結合距離の組み合わせに対して、高精度に原子パラメータ値を仮決定できるようになると考えられる。したがって、分子数が大きい場合には、パラメータ決定過程の終盤で選択された分子に関しては、その分子を構成するすべての原子に対して精度の高い原子パラメータ値が仮決定する確率が高くなる。これは、数９における分子nに関する和の順序を、原子パラメータ決定に用いる分子の計算順序と対応させた場合に、nが小さい時にはδ(n)が１となる確率が高く、nが大きくなるにしたがってδ(n)が0となる確率が高くなることに相当する。分子数Nが小さければ、パラメータ決定過程の終盤であってもδ(n)が０となる確率はそれほど高くはないと考えられるが、分子数Nが大きくなれば、δ(n)が０となる確率が高まると考えられる。
【００８１】
また、本実施例の方法は、大きな分子の計算回数を削減する効果もある。これまでは、原子パラメータ決定に要する計算量を示す数２と数９の比較において、f(A(n))に掛かる定数項部分の大小を議論してきた。しかし、f(A(n))は分子nの原子数A(n)のべき乗関数であるので、A(n)の大きな分子に関する計算回数の削減がより重要である。本実施例における分子の選択順序は、一般的に言うと、小さな分子が先、大きな分子が後となっている。最初に選択する分子は、原子パラメータの数が１である分子であり、これにはメタン、アンモニアなど非水素原子数が１の分子が相当する。ベンゼン分子のように対称性が高いために６つの炭素原子を含むにも関わらず、原子パラメータの数が１となる例外もあるが、通常は、小さな分子が優先的に選択され、大きな分子はパラメータ決定過程の終盤で選択される。先ほど議論したように、パラメータ決定過程の終盤では、δ(n)が０となる確率が高い。つまり、未決定の原子パラメータを決定するための計算の反復を必要とする確率が高い。
【００８２】
以上の議論を図２６にまとめる。原子パラメータ決定に用いる分子数Nが十分に大きく、また、分子は数２、数９の和記号に対して本実施例における計算順に並んでいるとする。計算の序盤2601、すなわちnが小さい時と、計算の終盤2602、すなわちnが大きいときに関わらず、従来法による分子nの計算量2603は一定である。全ての分子の溶媒和エネルギーおよび原子パラメータ値の変化に対する溶媒和エネルギー勾配計算を毎回反復する手法であるため、分子の計算順序の変更は計算量に影響を与えない。一方で、本実施例による分子nの計算量2604は、計算の序盤では、未決定パラメータの決定処理を実行する可能性が高いため、式８においてδ(n)=1とした値になる。しかし、計算の序盤2601では選択される分子は小さい可能性が高いため、結果として一回の溶媒和エネルギー計算に要する計算量f(A(n))2605は小さい。計算の終盤における分子nの計算量2604は、分子を構成する全ての原子に対して原子パラメータが仮決定している可能性が高いので、数９においてδ(n)=0とした値になる。計算の終盤2601では選択される分子は大きい可能性が高いため、結果として一回の溶媒和エネルギー計算に要する計算量f(A(n))2605は大きい。
（実施例５）
本実施例では、６３分子を用いた許容誤差固定による原子パラメータの決定結果を示す。原子パラメータ決定に用いた分子の内訳は、陽イオン１８個、陰イオン１１個、中性分子３４個である。許容エネルギー誤差は、陽イオン、陰イオンについては、1.4 kcal/mol、中性分子については0.20 kcal/molと設定した。結果として、５７個の原子パラメータを得た。
【００８３】
図２７に本実施例と従来法における計算量の比較を示す。本実施例では、５７個の原子パラメータを逐次決定する過程において、各分子の溶媒和エネルギー計算の反復回数2701は平均して3.32回であった。また、一台のコンピュータを用いた原子パラメータ決定に要した計算時間で定義する計算量2702は約504分であった。勾配行列を用いたときの計算量は実測できなかったが、以下の計算により計算量の下限値2703を3183分と見積もった。まず、数６において、各分子nについての一回の溶媒和エネルギー計算に必要な計算量f(A(n))には実測定値が存在する。また、各分子に含まれる原子タイプの数は、Tomasiらによる原子タイプ定義[V. Barone, M. Cossi, and J. Tomasi, Journal of Chemical Physics, Vol 107, 3210, (1997)] を用いて決定した。残された未知の値は反復回数Iのみである。行列方程式により定義する目的関数の最小化は多変量非線形問題であり、求解に要する反復回数を削減する方法は数多くある。しかし、反復回数Iは、一変数方程式の求解に要する反復回数の平均値である3.32回よりは多いと考えられる。そこで、反復回数を3.32以上の最小の整数I=4とし、見積もりを得た。この結果から、本実施例による原子パラメータ決定方法は、従来法よりも最低でも約6.3倍計算量が少ないと言える。従来法では、所望の精度を得られない場合には、原子タイプを定義し直し、原子パラメータ決定の全過程を反復する必要があるので、この反復回数をkとすると、計算量の差は6.3k倍となる。
（実施例６）
本実施例の許容誤差固定によるパラメータ決定方法は、分子を「解析対象物」、分子を構成する要素を「解析対象物を構成する要素」、要素パラメータを「要素パラメータ」、溶媒和エネルギーの実験値を単に「実験値」、溶媒和エネルギーの計算値を単に「計算値」とし、一般化可能である。以下に、図２８、図２９を用いて、一般化された許容誤差固定によるパラメータ決定方法を説明する。
【００８４】
図２８は、許容誤差固定による要素パラメータ決定方法のフローチャートである。まず、入力装置2801から許容誤差ε28011、収束閾値τ28012、解析対象物の計算に必要な情報および実験値28013、解析対象物構成する各要素の要素タイプ、および、各要素の要素属性値28014を入力する。次に、初期化ステップ2802では、要素タイプごとに設定した要素パラメータの初期化ステップ28021と、全解析対象物を未決定要素パラメータ数の昇順に並べ替えるステップ28022を実行する。計算順序リスト生成ステップ2803では、未決定の要素パラメータを含む解析対象物から、未決定の要素パラメータ数が１以下である解析対象物を選択するステップ28031と、選択解析対象物の並び順を計算順序として決定し、計算順序リストを生成するステップ28032を実行する。次に、一変数方程式求解ステップ2804では、計算順序リストの最上位解析対象物が未決定要素パラメータを含むかの判定ステップ28041で含むと判定された場合、計算値の変化量がτ未満となるまで、一変数方程式の求解を反復し、未決定要素パラメータを決定するステップ28042を実行する。つぎに、判定ステップ28041の結果に関わらず、要素パラメータを用いて該解析対象物の計算値と実験値の差である計算誤差を算出するステップ28043と、計算誤差の絶対値と許容誤差を比較するステップ28404を実行し、該計算誤差の方が許容誤差を超える場合には、仮決定の要素パラメータの１つを未決定の要素パラメータとして新規に設定するステップ28045を経てステップ28041に戻り、計算誤差の方が許容誤差以下の場合には計算順序リスト更新ステップ28046を実行する。計算順序リスト終了判定ステップ28047で、計算順序リストに最上位の解析対象物が存在する場合にはステップ28041に戻る。存在しない場合には、原子パラメータ更新ステップ2805へと進む。次に、要素パラメータ更新ステップ2805では、決定した要素パラメータに対応する要素タイプに関して、要素属性値電の関数である要素パラメータの表現式を更新するステップ28051と、未決定の要素の一部または全部の要素パラメータ値を仮決定するステップ28052と、未決定解析対象物を未決定要素パラメータの昇順に並べ替えるステップ28053と、未決定要素パラメータ数が同一の解析対象物を要素パラメータ数の降順に並べ替えるステップ28054と、未決定解析対象物を更新するステップ28055を実行する。要素パラメータ計算終了判定ステップ2806では、未決定解析対象物が存在するか判定する。存在する場合には、計算順序リスト生成ステップ2803、一変数方程式求解ステップ2804、要素パラメータ更新ステップ2805を繰り返す。存在しない場合には、出力装置2807より、決定した要素パラメータ値28071を出力する。
【００８５】
図２９は、要素パラメータ決定方法における、計算順序リスト生成ステップ、一変数方程式求解ステップ、および要素パラメータ更新ステップにおける処理の詳細である。例として、解析対象物数が９、要素タイプ数が８の場合における各ステップの一回目と二回目の繰り返し過程を説明する。計算順序リスト生成ステップ2901で実行する一回目の処理3202では、あらかじめ未決定要素パラメータ数29021の昇順に並べられている解析対象物から、未決定要素パラメータ数が１以下である解析対象物29022（解析対象物８、解析対象物３、解析対象物４）を選択する。これらの解析対象物の計算順序は、解析対象物の並び順と同じ解析対象物８、解析対象物３、解析対象物４の順であり、この順序を計算順序リスト29023に格納する。一変数方程式求解ステップ2903で実行する一回目の処理2904では、計算順序リストに基づき未決定要素パラメータを決定する。最初の計算は解析対象物８であるので、要素パラメータra,8を用いて計算値を取得し、実験値E8を再現できるようにra,8を決定する。この過程では、ra,8を変化させつつ実験値を再現できる値を探索するので、一変数方程式の反復求解が必要となる。選択された残り二つの解析対象物についても同様の処理を実行する。要素パラメータ更新ステップ2905で実行する一回目の処理2906では、決定した要素パラメータが対応する要素タイプ29061であるa, b, およびhに関して、要素属性値の関数である要素パラメータの表現式を更新する。要素パラメータの表現式は、要素属性値を変数に持たない関数であってもよい。例えば、要素タイプaの要素パラメータ29062はra(p1,p2)=c0=ra,8と表現できる。要素タイプaに分類されている要素の要素パラメータ値は、ra(p1,p2)=c0と仮決定できる。ただし、要素属性値がパラメータ決定に用いられた要素の要素属性値とある閾値以上離れた要素29063に関しては、要素パラメータ値を仮決定しない選択肢もある。要素パラメータの値が決定している要素2910を●、仮決定している要素2911を▲、値が割り当てられていない、すなわち未決定要素2912を○で示す。次に、この時点で未決定の解析対象物である、解析対象物５、解析対象物６、解析対象物７、解析対象物２、解析対象物１、解析対象物９の６つを未決定要素パラメータ数29064の昇順に並べ替える。また、未決定要素パラメータ数が等しい場合には、要素パラメータ数29065の降順に並べ替える。これら６つの解析対象物全てが、新たな未決定解析対象物となる。計算順序リスト生成ステップ2901で実行する二回目の処理2907では、あらかじめ並べられている解析対象物から、未決定要素パラメータ数29071が１以下である解析対象物29072（解析対象物２、解析対象物５、解析対象物６）を選択する。これらの解析対象物の計算順序は解析対象物２、解析対象物５、解析対象物６の順であり、この順序を計算順序リスト29073に格納する。一変数方程式求解ステップ2903で実行する二回目の処理2908では、計算順序リストに基づいた未決定要素パラメータを決定する。未決定パラメータが存在しない解析対象物２の場合について特に説明する。まず、与えられた要素パラメータを用いて計算値を取得する。次に、計算値と実験値の差と定義する計算誤差の絶対値が許容誤差未満となっているかを確認する。計算誤差の絶対値が許容誤差以上の場合は、解析対象物２に含まれる二つ要素パラメータra,2とrh,2のうちいずれかを未決定要素パラメータとして新規に設定する。選択する要素パラメータは、両要素パラメータの変化に対する計算値の勾配が大きい方とする。選択されたパラメータがra,2である場合には、実験値E2を再現するようにra,2の値を決定する。要素パラメータ更新ステップ2905で実行する二回目の処理2909では、決定した要素パラメータが対応する要素タイプ29091であるa, b, h, dおよびeに関する、要素属性値の関数である要素パラメータの表現式を更新する。例えば、要素タイプaの要素パラメータ29092はra(p1,p2)=c0=ra,8であったが、解析対象物２の計算結果から、要素タイプaの要素パラメータには新たに決定したra,2の値も用いる必要があるので、より複雑な関数形を用いる必要がある。例えば、ra(p1,p2)=c'0+c'1p1と更新してもよい。表現式を更新した後に、該当する要素タイプに分類される要素の要素パラメータ値を決定し、未決定の解析対象物を並べ替える操作を反復の一回目同様に実行する。この処理を未決定の解析対象物がなくなるまで反復する。
【００８６】
最後に従来法による要素パラメータの決定方法と、本実施例による要素パラメータ決定方法における計算量の比較を行う。従来法における計算量は式１に記載している。本実施例では、解析対象物を順次選択し、各解析対象物について計算値の変化量が収束閾値未満となるまで計算を反復する。解析対象物ごとに異なる反復回数をI(n)とする。但し、未決定要素パラメータ数が０である解析対象物も存在する。そこで、未決定要素パラメータ決定処理を必要とする解析対象物では１、必要としない解析対象物では０と定義するδ関数を用いると、計算量は、（数１０）
【００８７】
【数１０】
2007323099.tif 000012
となる。
【００８８】
数２と数１０を比較すると、二つの点で、本実施例の方法は従来法よりも計算量が少ないことがわかる。第一に、従来法では式１の２行目右辺における括弧内の第一項が”I”であるのに対し、本実施例では括弧内の第一項が”1”であるため計算量が少ない。第二に、従来法では括弧内の第二項が”I・p(n)”であるのに対し、本実施例では括弧内の第一項が”I(n)・δ(n)”であるため以下の二つの理由により、計算量が少ない。一つ目に、多変量非線形問題である行列方程式の求解に要する反復回数Iは、一変量非線形問題である一変数方程式の求解に要する反復回数I(n)よりも大きいと考えられる。二つ目に、解析対象物nに含まれる要素タイプ数p(n)の最低値が１であるのに対し、δ関数は０か１の値である。
（実施例７）
図３０〜図３３を用いて、原子パラメータ決定方法に付随するユーザーインターフェースの実装例を説明する。
【００８９】
図３０は、本実施例の原子パラメータ決定方法における、情報表示ステップを含むフローチャートである。情報表示に関連しない箇所については、実施例１、実施例４において詳細に説明したため、ここでは簡略化して示している。本実施例による原子パラメータ決定方法では、まず、入力ステップ3001と初期化ステップ3002を実行する。次に、計算順序リスト生成ステップ3003では、内部処理30031を経て、計算順序リスト表示ステップ30032で、計算順序リストを画面表示する。ユーザーは画面表示された計算順序に変更を加えてもよく、最終的に計算順序を確定する。計算進行状況表示ステップ30033では、原子パラメータ計算の進行状況の概要と、分子ごとの計算進行状況を表示する。次に、一変数方程式求解ステップ3004を実行する。原子パラメータ更新ステップ3005では、内部処理30051を経て、原子パラメータ表示ステップ30052では、原子パラメータを、原子パラメータ、原子電荷、および平均結合距離の三者が形成する空間を用いて表示する。最後に、終了判定ステップ3006を経て、出力ステップ3007を実行する。以上のフローチャート中で、情報表示に関係するのは、計算順序リスト表示ステップ30032、計算進行状況表示ステップ30033、原子パラメータ表示ステップ30052の三箇所である。
【００９０】
図３１は、計算順序リスト表示ステップにおける画面表示例である。表示対象物表示タブ3101に含まれる三つの表示対象物のうち、計算順序リストタブ31011が選択されている。計算順序リスト表示画面は、計算順序表示領域3102、分布図表示領域3103、分子構造図表示領域3104に三分割されている。表示の初期段階では、計算順序表示領域3102に、計算順序リスト31021が表示されている。表は行ごとに選択可能であり、図中においては分子５が選択状態になっている。選択した分子の計算順序は、上矢印ボタン31022や下矢印ボタン31023により調節できる。また、削除ボタン31024により、選択分子の計算を行わないという選択もできる。ユーザーが確定ボタン31025を押すと、計算順序が確定する。ユーザーが計算順序の変更を判断するための補助データが、分布図と、分子構造図である。分布図表示領域3103には、分子５の未決定原子パラメータが属する原子タイプdに分類される原子の分布図の例を示す。原子は、原子電荷と平均結合距離が形成する平面上にプロットする。プロットは、各原子に対応する原子パラメータの「値あり」、「値なし」などの状態や、計算順序リスト19021における選択状態を反映して、異なる形状あるいは模様を有する。本例では、原子タイプdに分類される原子は二つしかなく、両者は、平面上の離れた位置にプロットされている。したがって、一方の原子パラメータを決定しても、決定した原子パラメータに基づきもう一方の原子パラメータ値を高精度に仮決定できるとは限らないことが視覚的に分かる。この情報に基づき、ユーザーが、分子５の計算の優先順位を下げる選択をしてもよい。分子構造図表示領域3104には、分子の立体構造と、分子を形成する各原子の原子パラメータ値に対応する半径を持つ球を重ねて表示する。ここでも、各原子に対応する原子パラメータの「値あり」、「値なし」などの状態、および計算順序リスト31021における選択状態を反映して、球を異なる模様で表示する。分子５は三つの非水素原子からなる分子であり、そのうち二つは対称性の観点から等価であると分かる。原子Xに対する原子パラメータ値は仮決定しており、値が割り当てられている。一方、原子Yの原子パラメータは、これから決定することが視覚的に分かる。分布図表示領域3103中のプロットや、分子構造図表示領域3104に表示された分子中の原子も選択可能である。プロットを選択すると、対応する原子が属する分子の構造を分子構造図表示領域3104に表示する。また、原子を選択すると、対応する原子タイプに関する原子の分布図を分布図表示領域3103に表示する。これらの機能により、着目している分子５のパラメータ決定に直接関連する他の分子の情報を取得可能となる。
【００９１】
図３２は、計算進行状況表示ステップにおける画面表示例である。表示対象物表示タブ3201に含まれる三つの表示対象物のうち進行状況タブ32011が選択されている。計算進行状況表示画面は、原子パラメータ計算進行状況の概要表示領域3202と各分子の原子パラメータ計算進行状況表示領域3203に分割されている。原子パラメータ計算進行状況の概要表示領域3202には、「決定」、「計算中」、「未決定」といった分子の計算状態の集計結果、および計算の経過時間を表示する。また、各分子の原子パラメータ計算進行状況表示領域3203には、各分子の計算状態、分子ID、原子パラメータ数、未決定原子パラメータ数を表示する。これらの情報から、パラメータ決定手順の進行度合いを把握可能となる。
【００９２】
図３３は、原子パラメータ表示ステップにおける画面表示例である。表示対象物表示タブ3301に含まれる三つの表示対象物のうち原子パラメータタブ33011が選択されている。原子パラメータ表示画面は、各原子タイプの原子パラメータ決定状況表示領域3302と、原子パラメータ値表示領域3303に分割されている。表示の初期段階では、各原子タイプの原子パラメータ決定状況表示領域3302に、原子タイプ、該当原子タイプに属する全原子数、原子パラメータが決定した原子数、仮決定の原子数、未決定の原子数を表形式で表示する。表は行ごとに選択可能であり、図中においては原子タイプaが選択状態になっている。原子パラメータ値表示領域3303には、選択された原子タイプにおける原子パラメータ値を、原子パラメータ、原子電荷、平均結合距離が作る空間で表示する。原子パラメータ値の三次元空間におけるグラフ表示により、原子パラメータの原子電荷と平均結合距離に対する依存性を直感的に理解できる。
【００９３】
以下、これまでに記載の実施例の効果について記載する。
【００９４】
本発明の実施例による、高確度に決定できる原子パラメータを含む分子を優先的に選択し、選択された分子に関する溶媒和エネルギー計算の順序を計算順序リストにより制御する方式により、従来よりも少ない計算量で、原子パラメータ決定が可能となる。
【００９５】
実施例１では、分子の溶媒和エネルギー計算に用いる原子パラメータを従来法よりも少ない計算量で決定可能であることを示した。実施例２では、１０分子を用いた実証計算により、本実施例の原子パラメータ決定方法は、従来法に基づく方法よりも少ない計算量で同等の精度を持つ原子パラメータを決定可能であることを示した。実施例４では、許容誤差固定の原子パラメータ決定方法により、分子の溶媒和エネルギー計算に用いる原子パラメータを従来法よりも少ない計算量で決定でき、かつ計算誤差も制御できることを示した。実施例５では、６３分子を用いた実証計算により、許容誤差固定による原子パラメータ決定方法に要する計算量は、従来法に基づく方法に要する計算量の見積もりよりも少ないことを示した。
【００９６】
また本発明の実施例による、高確度に決定できるパラメータを含む解析対象物を優先的に選択し、選択された解析対象物に関する計算を計算順序リストにより制御する方式により、従来よりも少ない計算量で、パラメータ決定が可能となる。
【００９７】
本発明の実施例は、計算対象が要素に分割でき、各要素に属性値が定義でき、かつ、計算値や計算値のパラメータ変化に対する勾配の評価が計算の律速段階である場合に適用可能である。
【００９８】
さらに、本発明の実施例では、パラメータ決定方法の重要な要素である解析対象物の計算順序の決定過程を可視化し、ユーザーが計算順序に変更を加える機能を設けたため、ユーザーの知識をパラメータ決定に反映可能である。実施例７では、ユーザーが計算順序の変更を判断するための補助データの表示例を示した。
【図面の簡単な説明】
【００９９】
【図１】連続誘電体モデルと、そこに導入するパラメータ。
【図２】勾配行列を用いた原子パラメータ決定方法のフローチャート。
【図３】勾配行列を用いた原子パラメータ決定方法における、行列生成ステップ、行列方程式求解ステップ、および原子パラメータ更新ステップにおける処理の詳細。
【図４】勾配行列を用いた原子パラメータ決定方法における原子タイプ定義、得られた原子パラメータ値、および各分子に対応する原子タイプ。
【図５】勾配行列を用いた原子パラメータ決定方法における計算量および計算誤差。
【図６】勾配行列を用いたパラメータ決定方法のフローチャート。
【図７】勾配行列を用いたパラメータ決定方法における、行列生成ステップ、行列方程式求解ステップ、および要素パラメータ更新ステップにおける処理の詳細。
【図８】原子パラメータ決定方法のフローチャート。
【図９】原子パラメータ決定方法における、計算順序リスト生成ステップ、一変数方程式求解ステップ、および原子パラメータ更新ステップにおける処理の詳細。
【図１０】分子に帰属する情報の入力。
【図１１】原子に帰属する情報の入力。
【図１２】原子パラメータの初期化。
【図１３】一変数方程式求解ステップにおけるハードウェアの制御方法。
【図１４】原子パラメータ値を、一つの原子属性値を変数とする一次関数で表現する場合の原子パラメータの更新例 。
【図１５】原子パラメータ値を、二つの原子属性値を変数とする一次関数で表現する場合の原子パラメータの更新例。
【図１６−１】本発明の実施例の適用に用いた原子タイプ定義、原子パラメータ値の決定過程、および各分子に対応する原子タイプ。
【図１６−２】本発明の実施例の適用に用いた原子タイプ定義、原子パラメータ値の決定過程、および各分子に対応する原子タイプ。
【図１７】従来法と本発明の実施例の計算量および計算誤差の比較。
【図１８】溶媒和エネルギー予測精度の比較。
【図１９】パラメータ決定方法のフローチャート。
【図２０】パラメータ決定方法における、計算順序リスト生成ステップ、一変数方程式求解ステップ、および要素パラメータ更新ステップにおける処理の詳細。
【図２１】許容誤差固定による原子パラメータ決定方法のフローチャート 。
【図２２】許容誤差固定による原子パラメータ決定方法における、計算順序リスト生成ステップ、一変数方程式求解ステップ、および原子パラメータ更新ステップにおける処理の詳細 。
【図２３】計算順序リスト更新ステップの詳細。
【図２４】原子パラメータ値を、一つの原子属性値を変数とする関数で表現する場合の原子パラメータの更新例 。
【図２５】原子パラメータ値を、二つの原子属性値を変数とする関数で表現する場合の原子パラメータの更新例 。
【図２６】原子パラメータ決定に用いる分子数Nが十分に大きいときの、計算量の比較。
【図２７】６３分子を用いた原子パラメータ決定における、本発明の実施例と従来法における計算量の比較 。
【図２８】許容誤差固定によるパラメータ決定方法のフローチャート 。
【図２９】許容誤差固定によるパラメータ決定方法における、計算順序リスト生成ステップ、一変数方程式求解ステップ、および要素パラメータ更新ステップにおける処理の詳細 。
【図３０】原子パラメータ決定方法における情報表示ステップを含むフローチャート。
【図３１】計算順序リスト表示ステップにおける画面表示例。
【図３２】計算進行状況表示ステップにおける画面表示例。
【図３３】原子パラメータ表示ステップステップにおける画面表示例。
【符号の説明】
【０１００】
101            計算対象分子
102            水溶媒
103            分子が占有する空間
104            誘電体が占有する空間
105            原子上に配置された球
106            原子上に配置された球の半径
201            入力装置
2011          収束エネルギー閾値τの入力
2012          分子の構造および溶媒和エネルギーの実験値の入力
2013          分子を構成する各原子の原子タイプの入力
202            原子パラメータ初期化ステップ
203            行列生成ステップ
204            行列方程式求解ステップ
205            原子パラメータ更新ステップ
206            原子パラメータ計算終了判定ステップ
207            出力装置
2071          収束した原子パラメータ値の出力
301            行列生成ステップ
302            行列生成ステップで実行する処理
3021          全分子の溶媒和エネルギー
3022          各原子パラメータの変化に対する全分子の溶媒和エネルギー勾配
3023          行数が分子数、列数が（原子タイプ数＋１）からなる行列
3024          行列要素
303            行列方程式求解ステップ
304            行列方程式求解ステップで実行する処理
3041          行列
3042          （原子タイプ数＋１）個の要素数を持つ解ベクトル
3043          分子数と同一の要素数を持つ溶媒和エネルギー実験値ベクトル
3044          目的関数
305            原子パラメータ更新ステップ
306            原子パラメータ更新ステップで実行する処理
3061          新しい原子パラメータベクトル
3062          更新前の原子パラメータベクトル
3063          行列方程式求解ステップ303で得られた解ベクトルから先頭要素のみを削除したベクトル
401           原子タイプ
402           得られた原子パラメータ値
403            各原子が含有する原子タイプ
501            計算量
502            計算誤差
503            ひし形プロット
601            入力装置
6011          収束閾値τの入力
6012          解析対象物の計算に必要な情報および実験値の入力
6013          解析対象物を構成する各要素の要素タイプの入力
602            要素パラメータ初期化ステップ
603            行列生成ステップ
604            行列方程式求解ステップ
605            要素パラメータ更新ステップ
606            要素パラメータ計算終了判定ステップ
607            出力装置
6071          収束した要素パラメータ値の出力
701            行列生成ステップ
702            行列生成ステップで実行する処理
7021          全解析対象物の計算値
7022          各要素パラメータの変化に対する全解析対象物の計算値勾配
7023          行数が解析対象物数、列数が（要素タイプ数＋１）からなる行列
7024          行列要素
703            行列方程式求解ステップ
704            行列方程式求解ステップで実行する処理
7041          行列
7042          （要素タイプ数＋１）個の要素数を持つ解ベクトル
7043          解析対象物数と同一の要素数を持つ実験値ベクトル
7044          目的関数
705            要素パラメータ更新ステップ
706            要素パラメータ更新ステップで実行する処理
7061          新しい要素パラメータベクトル
7062          更新前の要素パラメータベクトル
7063          行列方程式求解ステップ303で得られた解ベクトルから先頭要素のみを削除したベクトル
801            入力装置
8011          収束エネルギー閾値τの入力
8012          分子の三次元座標情報および溶媒和エネルギーの実験値の入力
8013          各原子の原子タイプ、および、各原子の原子電荷と平均結合距離の入力
802            初期化ステップ
8021          原子パラメータの初期化ステップ
8022          全分子の未決定原子パラメータ数の昇順への並べ替えステップ
803            計算順序リスト生成ステップ
8031          分子選択ステップ
8032          計算順序決定ステップ
804            一変数方程式求解ステップ
8041          未決定原子パラメータ決定ステップ
8042          計算順序リスト更新ステップ
8043          計算順序リスト終了判定ステップ
805            原子パラメータ更新ステップ
8051          原子パラメータ表現式更新ステップ
8052          未決定分子の未決定原子パラメータの昇順への並べ替えステップ
8053          未決定原子パラメータ数同一分子の原子パラメータ数の降順への並べ替えステップ
8054          未決定分子を更新するステップ
806            原子パラメータ計算終了判定ステップ
807            計算誤差算出ステップ
808            出力装置
8081          決定した原子パラメータ値の出力
901            計算順序リスト生成ステップ
902            計算順序リスト生成ステップで実行する一回目の処理
9021          未決定原子パラメータ数
9022          未決定原子パラメータ数が１である分子
9023          計算順序リスト
903            一変数方程式求解ステップ
904            一変数方程式求解ステップで実行する一回目の処理
905            原子パラメータ更新ステップ
906            原子パラメータ更新ステップで実行する一回目の処理
9061          決定した原子パラメータが対応する原子タイプ
9062          原子タイプaの原子パラメータ
9063          原子電荷が同一原子タイプに属する他の原子の原子電荷とある閾値以上離れた原子
907            計算順序リスト生成ステップで実行する二回目の処理
9071          未決定原子パラメータ数
9072          未決定原子パラメータ数が１である分子
9073          計算順序リスト
908            一変数方程式求解ステップで実行する二回目の処理
909            原子パラメータ更新ステップで実行する二回目の処理
9091          決定した原子パラメータが対応する原子タイプ
9092          原子タイプaの原子パラメータ
910            「値あり」の原子
911            「値なし」（未決定）の原子
1001          分子の三次元座標情報
1002          パラキシレン分子の構造
1003          溶媒和エネルギー実験値
1004          パラキシレン分子の溶媒和エネルギーの実験値
1005          原子パラメータ数
1006          パラキシレン分子の原子パラメータ数
1101          原子タイプ
1102          原子電荷
1103          平均結合距離
1201          原子に帰属する情報
1202          原子パラメータの値
1203          原子パラメータの状態
1204          原子パラメータの状態（値あり）
1205          原子パラメータの状態（値なし）
1206          分子に帰属する情報
1207          分子に含まれる「値あり」の原子パラメータ数
1208          分子に含まれる「値なし」の原子パラメータ数
1209          分子の計算状態
1210          分子に帰属する情報の初期値
1211          分子に含まれる「値あり」の原子パラメータ数の初期値
1212          分子に含まれる「値なし」の原子パラメータ数の初期値
1213          分子の計算状態の初期値
1301          計算順序リスト
1302          入力ファイル生成ステップ
1303          溶媒和エネルギー計算用の入力ファイル
1304          溶媒和エネルギー計算ジョブ投入ステップ
1305          計算機
1306          溶媒和エネルギー計算プログラム
1307          計算結果出力ステップ
1308          溶媒和エネルギー計算結果の出力ファイル  
1309          出力結果処理ステップ
1401          原子パラメータ値
1402          原子属性値
1403          原子パラメータが「値あり」の原子
1404          原子パラメータが「値なし」の原子
1405          原子パラメータ値と原子属性値を関連付ける一次関数
1406          「値なし」の原子の原子属性値
1407          「値あり」の原子の原子属性値と大きくかけ離れた原子属性値を持つ「値なし」の原子
1501          原子パラメータ値
1502          原子電荷
1503          平均結合距離
1504          原子パラメータが「値あり」の原子
1505          原子パラメータが「値なし」の原子
1506          原子パラメータ値と原子属性値を関連付ける一次関数
1507          「値なし」の原子の原子電荷と平均結合距離
1508          「値あり」の原子の原子属性値と大きくかけ離れた原子属性値を持つ「値なし」の原子
1601          原子タイプ
1602          原子パラメータ値の決定過程
1603          各原子が含有する原子タイプ
1604          未決定原子パラメータ数が1である分子
1605          原子パラメータrCH3
1606          原子パラメータrCH3
1607          未決定の原子パラメータ
1608          未決定原子パラメータ数
1609          未決定原子パラメータ数が1である分子
1610          原子タイプNH2に関して決定した４つの原子パラメータ
1611          原子パラメータrNH2
1701          計算量
1702          計算誤差
1703          ひし形プロット
1704          長方形プロット
1801          分子
1802          溶媒和エネルギーの実験値
1803          従来法による計算誤差
1804          本発明の実施例による計算誤差
1805          平均絶対誤差
1901          入力装置
19011        収束エネルギー閾値τの入力
19012        解析対象物の計算に必要な情報および実験値の入力
19013        各要素の要素タイプ、および、各要素の要素属性値の入力
1902          初期化ステップ
19021        要素パラメータの初期化ステップ
19022        全解析対象物の未決定要素パラメータ数の昇順への並べ替えステップ
1903          計算順序リスト生成ステップ
19031        解析対象物選択ステップ
19032        計算順序決定ステップ
1904          一変数方程式求解ステップ
19041        未決定要素パラメータ決定ステップ
19042        計算順序リスト更新ステップ
19043        計算順序リスト終了判定ステップ
1905          要素パラメータ更新ステップ
19051        要素パラメータ表現式更新ステップ
19052        未決定解析対象物の未決定要素パラメータの昇順への並べ替えステップ
19053        未決定要素パラメータ数同一解析対象物の要素パラメータ数の降順への並べ替えステップ
19054        未決定解析対象物を更新するステップ
1906          要素パラメータ計算終了判定ステップ
1907          計算誤差算出ステップ
1908          出力装置
19081        決定した要素パラメータ値の出力
2001          計算順序リスト生成ステップ
2002          計算順序リスト生成ステップで実行する一回目の処理
20021        未決定要素パラメータ数
20022        未決定要素パラメータ数が１である解析対象物
20023        計算順序リスト
2003          一変数方程式求解ステップ
2004          一変数方程式求解ステップで実行する一回目の処理
2005          要素パラメータ更新ステップ
2006          要素パラメータ更新ステップで実行する一回目の処理
20061        決定した要素パラメータが対応する要素タイプ
20062        要素タイプaの要素パラメータ
20063        要素電荷が同一要素タイプに属する他の要素の要素電荷とある閾値以上離れた要素
2007          計算順序リスト生成ステップで実行する二回目の処理
20071        未決定要素パラメータ数
20072        未決定要素パラメータ数が１である解析対象物
20073        計算順序リスト
2008          一変数方程式求解ステップで実行する二回目の処理
2009          要素パラメータ更新ステップで実行する二回目の処理
20091        決定した要素パラメータが対応する要素タイプ
20092        要素タイプaの要素パラメータ
2010          「値あり」の要素
2011          「値なし」（未決定）の要素
2101          入力装置
21011        許容エネルギー誤差εの入力
21012        収束エネルギー閾値τの入力
21013        分子の三次元座標情報および溶媒和エネルギーの実験値の入力
21014        各原子の原子タイプ、および、各原子の原子電荷と平均結合距離の入力
2102          初期化ステップ
21021        原子パラメータの初期化ステップ
21022        全分子の未決定原子パラメータ数の昇順への並べ替えステップ
2103          計算順序リスト生成ステップ
21031        分子選択ステップ
21032        計算順序決定ステップ
2104          一変数方程式求解ステップ
21041        計算順序リストの最上位分子が未決定原子パラメータを含むかの判定ステップ
21042        未決定原子パラメータ決定ステップ
21043        計算誤差算出ステップ
21044        計算誤差と許容エネルギー誤差の比較ステップ
21045        未決定原子パラメータ新規設定ステップ
21046        計算順序リスト更新ステップ
21047        計算順序リスト終了判定ステップ
2105          原子パラメータ更新ステップ
21051        原子パラメータ表現式更新ステップ
21053        未決定分子の未決定原子パラメータの昇順への並べ替えステップ
21054        未決定原子パラメータ数同一分子の原子パラメータ数の降順への並べ替えステップ
21055        未決定分子を更新するステップ
2106          原子パラメータ計算終了判定ステップ
2107          出力装置
21071        決定した原子パラメータ値の出力
2201          計算順序リスト生成ステップ
2202          計算順序リスト生成ステップで実行する一回目の処理
22021        未決定原子パラメータ数
22022        未決定原子パラメータ数が１以下である分子
22023        計算順序リスト
2203          一変数方程式求解ステップ
2204          一変数方程式求解ステップで実行する一回目の処理
2205          原子パラメータ更新ステップ
2206          原子パラメータ更新ステップで実行する一回目の処理
22061        決定した原子パラメータが対応する原子タイプ
22062        原子タイプaの原子パラメータ
22063        原子電荷が同一原子タイプに属する他の原子の原子電荷とある閾値以上離れた原子
2207          計算順序リスト生成ステップで実行する二回目の処理
22071        未決定原子パラメータ数
22072        未決定原子パラメータ数が１以下である分子
22073        計算順序リスト
2208          一変数方程式求解ステップで実行する二回目の処理
2209          原子パラメータ更新ステップで実行する二回目の処理
22091        決定した原子パラメータが対応する原子タイプ
22092        原子タイプaの原子パラメータ
2210          「値あり」（決定）の原子
2211          「値あり」（仮決定）の原子
2212          「値なし」（未決定）の原子
2301          開始
2302          未決定パラメータの新規設定判定ステップ
2303          計算順序リストの全削除処理
2304          計算順序リストから最上位分子削除処理
2305          計算順序リスト終了判定ステップ
2306          終了（原子パラメータ更新ステップへ）
2307          終了（未決定原子パラメータ決定ステップへ）
2401          原子パラメータ値
2402          一つの原子属性値
2403          「値あり」の原子
2404          「値なし」の原子
2405          原子パラメータ値と原子属性値を関連付ける一次関数
2406          「値なし」の原子の原子属性値
2407          仮決定した原子パラメータ値がそのまま決定値となった原子
2408          仮決定した原子パラメータ値が未決定原子パラメータとして設定され、原子パラメータが再決定した原子
2409          原子属性値に基づいた三つの区域
2410          それぞれの区域に定義した一次関数
2501          原子パラメータ値
2502          原子電荷
2503          平均結合距離
2504          「値あり」の原子
2505          「値なし」の原子
2506          原子パラメータ値と原子属性値を関連付ける一次関数
2507          「値なし」の原子の原子電荷と平均結合距離
2508          仮決定した原子パラメータ値がそのまま決定値となった原子
2509          仮決定した原子パラメータ値が未決定原子パラメータとして設定され、原子パラメータが再決定した原子
2510          原子電荷と平均結合距離が成す平面上の二つの区域
2511          それぞれの区域に定義した一次関数
2601          計算の序盤
2602          計算の終盤
2603          従来法による分子nの計算量
2604          本発明の実施例による分子nの計算量
2605          一回の溶媒和エネルギー計算に要する計算量
2701          反復回数
2702          本発明の実施例を用いたときの計算量
2703          勾配行列を用いたときの計算量の下限値
2801          入力装置
28011        許容誤差εの入力
28012        収束閾値τの入力
28013        解析対象物の計算に必要な情報および実験値の入力
28014        各要素の要素タイプ、および、各要素の要素属性値電荷の入力
2802          初期化ステップ
28021        要素パラメータの初期化ステップ
28022        全解析対象物の未決定要素パラメータ数の昇順への並べ替えステップ
2803          計算順序リスト生成ステップ
28031        解析対象物選択ステップ
28032        計算順序決定ステップ
2804          一変数方程式求解ステップ
28041        計算順序リストの最上位解析対象物が未決定要素パラメータを含むかの判定ステップ
28042        未決定要素パラメータ決定ステップ
28043        計算誤差算出ステップ
28044        計算誤差と許容誤差の比較ステップ
28045        未決定要素パラメータ新規設定ステップ
28046        計算順序リスト更新ステップ
28047        計算順序リスト終了判定ステップ
2805          要素パラメータ更新ステップ
28051        要素パラメータ表現式更新ステップ
28053        未決定解析対象物の未決定要素パラメータの昇順への並べ替えステップ
28054        未決定要素パラメータ数同一解析対象物の要素パラメータ数の降順への並べ替えステップ
28055        未決定解析対象物を更新するステップ
2806          要素パラメータ計算終了判定ステップ
2807          出力装置
28071        決定した要素パラメータ値の出力
2901          計算順序リスト生成ステップ
2902          計算順序リスト生成ステップで実行する一回目の処理
29021        未決定要素パラメータ数
29022        未決定要素パラメータ数が１以下である解析対象物
29023        計算順序リスト
2903          一変数方程式求解ステップ
2904          一変数方程式求解ステップで実行する一回目の処理
2905          要素パラメータ更新ステップ
2906          要素パラメータ更新ステップで実行する一回目の処理
29061        決定した要素パラメータが対応する要素タイプ
29062        要素タイプaの要素パラメータ
29063        要素属性値が同一要素タイプに属する他の要素の要素属性値とある閾値以上離れた要素
2907          計算順序リスト生成ステップで実行する二回目の処理
29071        未決定要素パラメータ数
29072        未決定要素パラメータ数が１以下である解析対象物
29073        計算順序リスト
2908          一変数方程式求解ステップで実行する二回目の処理
2909          要素パラメータ更新ステップで実行する二回目の処理
29091        決定した要素パラメータが対応する要素タイプ
29092        要素タイプaの要素パラメータ
2910          「値あり」（決定）の要素
2911          「値あり」（仮決定）の要素
2912          「値なし」（未決定）の要素
3001          入力ステップ
3002          初期化ステップ
3003          計算順序リスト生成ステップ
30031        内部処理
30032        計算順序リスト表示ステップ
30033        計算進行状況表示ステップ
3004         方程式求解ステップ
3005          原子パラメータ更新ステップ
30051        内部処理
30052        原子パラメータ表示ステップ
3006          終了判定ステップ
3007          出力ステップ
3101          表示対象物表示タブ
31011        計算順序リストタブ
31021        計算順序リスト
31022        上矢印ボタン
31023        下矢印ボタン
31024        削除ボタン
31025        確定ボタン
3102          計算順序表示領域
3103          分布図表示領域
3104          分子構造図表示領域
3201         表示対象物表示タブ
32011        進行状況タブ
3202         原子パラメータ計算進行状況の概要表示領域
3203         各分子の原子パラメータ計算進行状況表示領域
3301 表示対象物表示タブ
33011        原子パラメータタブ
3302 原子パラメータ決定状況表示領域
3303 原子パラメータ値表示領域。
【図２】
2007323099.tif 000013
【図６】
2007323099.tif 000014
【図８】
2007323099.tif 000015
【図１１】
2007323099.tif 000016
【図１２】
2007323099.tif 000017
【図１３】
2007323099.tif 000018
【図１４】
2007323099.tif 000019
【図１８】
2007323099.tif 000020
【図１９】
2007323099.tif 000021
【図２１】
2007323099.tif 000022
【図２３】
2007323099.tif 000023
【図２４】
2007323099.tif 000024
【図２６】
2007323099.tif 000025
【図２７】
2007323099.tif 000026
【図２８】
2007323099.tif 000027
【図２９】
2007323099.tif 000028
【図３０】
2007323099.tif 000029
【図１】
2007323099.tif 000030
【図３】
2007323099.tif 000031
【図４】
2007323099.tif 000032
【図５】
2007323099.tif 000033
【図７】
2007323099.tif 000034
【図９】
2007323099.tif 000035
【図１０】
2007323099.tif 000036
【図１５】
2007323099.tif 000037
【図１６−１】
2007323099.tif 000038
【図１６−２】
2007323099.tif 000039
【図１７】
2007323099.tif 000040
【図２０】
2007323099.tif 000041
【図２２】
2007323099.tif 000042
【図２５】
2007323099.tif 000043
【図３１】
2007323099.tif 000044
【図３２】
2007323099.tif 000045
【図３３】
2007323099.tif 000046