非定常加熱による多層積層材の熱物性同定方法
【課題】積層試料を構成する物質の全ての熱物性値が未知の場合にも、熱物性値を同定することができる非定常加熱による多層積層材の熱物性同定方法を提供する。
【解決手段】４つの熱物性値を独立変数とした評価関数ＪをＪ＝ＷSＪS＋ＷPＪPと定義する。なお、ＷS、ＷPは重み係数、ＪSはステップ加熱した試料裏面の昇温曲線観測データと解析解との２乗誤差の積分値、ＪＰはフラッシュ加熱した試料裏面の昇温曲線観測データと解析解との２乗誤差の積分値である。第１層目を銅、第２層目をアルミとし、真値で計算した解析解を観測データとして用い、第１層の熱伝導率λ1、第２層の熱伝導率λ2を未知パラメータとして変化させた場合、評価関数Ｊは、図に示すようになり、多峰性が認められるが、極値が観測され、昇温曲線に解析解曲線群を評価関数が最小になるようにフィッティングさせることにより熱物性値を決定することができる。
【解決手段】４つの熱物性値を独立変数とした評価関数ＪをＪ＝ＷSＪS＋ＷPＪPと定義する。なお、ＷS、ＷPは重み係数、ＪSはステップ加熱した試料裏面の昇温曲線観測データと解析解との２乗誤差の積分値、ＪＰはフラッシュ加熱した試料裏面の昇温曲線観測データと解析解との２乗誤差の積分値である。第１層目を銅、第２層目をアルミとし、真値で計算した解析解を観測データとして用い、第１層の熱伝導率λ1、第２層の熱伝導率λ2を未知パラメータとして変化させた場合、評価関数Ｊは、図に示すようになり、多峰性が認められるが、極値が観測され、昇温曲線に解析解曲線群を評価関数が最小になるようにフィッティングさせることにより熱物性値を決定することができる。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、試料の熱伝導率、熱拡散率などの熱物性の同定方法、特に、多層積層材の各層の熱伝導率、熱拡散率などの熱物性を同定する方法に関する。
【背景技術】
【０００２】
金属、セラミックス、ガラス、カーボン、プラスチックなど個体材料の熱的特性を精度よく測定する方法としてフラッシュ加熱法がある。このフラッシュ加熱法は、試料の表面にレーザ光を瞬間的（フラッシュ状）に照射して加熱し、一定距離離れた場所に伝わる熱を温度変化から測定するものである。例えば、直径１０ｍｍ、厚さ１〜３ｍｍ程度の小円板状試料の片面をレーザ光で瞬時に加熱して反対側の温度挙動から熱測定を行う装置があり、熱的特性としては、熱拡散率、及び比熱容量を測定することができ、さらに密度が既知であれば、熱伝導率を算出することができる。
【０００３】
一方、近年の新しい材料開発の動向の一つとして、複数の素材を組み合わせることにより、今までに無かった新しい機能を持たせた材料を開発することが行われており、その中でも特に高熱伝導性基板などが注目されている。しかし、この分野で必要とされる評価測定の部門で、出来上がった新機能材料の各々の層の熱物性を直接測定する方法はまだ確立されていない。
【０００４】
このような多層試料の熱物性を測定する場合、従来、定常法での熱伝導率測定により、多層試料全体の熱伝導率を測定することが行われている。この測定法は試料板の片面を加熱、反対面を冷却し、試料板厚み方向に温度勾配をつけ、その温度差と試料厚み及び試料を貫通する熱流束から熱伝導率を算出する方法である。
【０００５】
また、各層の熱物性を測定する試みとして、上記のフラッシュ加熱法を使用することが提案されている。すなわち、未知の一層を含む多層材料の表面側からレーザフラッシュを照射して得られる裏面側の温度応答と、多層材料の熱物性による理論温度応答との２乗誤差を評価関数（ポテンシャル）とし、この評価関数が最小値になるパラメータを求めることにより高精度でかつ効率的に多層材料中の一層における未知の熱拡散率、比熱等の熱定数を決定するようにしている（例えば、特許文献１参照）。
【先行技術文献】
【特許文献】
【０００６】
【特許文献１】特開平８−２６１９６７号公報
【発明の概要】
【発明が解決しようとする課題】
【０００７】
上記のように、熱拡散率及び比熱が未知の一層を含む多層材料の表面側からレーザフラッシュを照射して得られる裏面側の温度応答と、熱拡散率、比熱、及びビオー数を変数として含む理論温度応答とを比較することにより、多層積層材の熱物性を同定することが提案されているが、この方法では、多層材料の一層の熱拡散率及び比熱が未知の場合のみ、熱物性を同定できるだけであり、積層試料を構成する物質の全ての熱物性値が未知の場合には同定することは困難である。
【０００８】
以下、２層積層試料を構成する各試料をそれぞれの形状（厚み）を既知として、未知の熱物性値の組み合わせについてフラッシュ加熱法を用いて最小２乗法で各パラメータ、すなわち熱物性値を決定できるかについて説明する。
未知パラメータを独立変数とした評価関数（ポテンシャル）を以下のＪPで定義し、このポテンシャルが単一の極値を持つなら、滑降シンプレックス法で観測された昇温曲線に解析解曲線群をフィッティングさせることによりパラメータ、すなわち、熱物性値を決定することができる。
ＪP(λ1,α1,λ2,α2)=∫{ＹP(t)-ｙP(t,λ1,α1,λ2,α2)}2dt
ただし、ＹP(t)はフラッシュ加熱で得られた観測信号データ、ｙP(t,λ1,α1,λ2,α2)はフィッティングさせる解析解を表す。
【０００９】
上記の解析解ｙP (t,λ1,α1,λ2,α2)は、以下の数式１のθr（ｔ）で表わされる。なお、数式１において、λ、α、ｄをそれぞれ熱伝導率、熱拡散率、試料厚みとしたとき、ｂは熱浸透率でｂ＝λ／√α、τは特性時間でτ＝ｄ2／αである。また、添え字１、２は１層目、２層目の物性値を表す。
【００１０】
【数１】
【００１１】
ただし、ｓkはｓを複素数として下記数式２のｋ番目の解である。
【００１２】
【数２】
【００１３】
ここで、第１層目（試料１）を銅（λ＝３９８Ｗ／ｍＫ、α＝１１７ｍｍ2／ｓ、厚み１ｍｍ）、第２層目（試料２）をアルミ（λ＝２３７Ｗ／ｍＫ、α＝９６．８ｍｍ2／ｓ、厚み１ｍｍ）とし、また、昇温波形の最大値は、試料からの熱放散がないと考え（すなわちビオー数をゼロとした）、観測データＹP(t)から知ることができるとし、λ1＝３９８ｍＫ、α1＝１１７ｍｍ2／ｓ、λ2＝２３７Ｗ／ｍＫ、α2＝９６．８ｍｍ2／ｓのときの上記数式１で計算した解析解を観測データＹP(t)として用い、試料１の熱伝導率λ1、試料２の熱伝導率λ2を未知パラメータとし、λ1を３９８±５０Ｗ／ｍＫ、λ2を２３７±５０Ｗ／ｍＫの範囲で変化させた場合のポテンシャルＪPは、図１に示すようになり、極値の識別は困難である。なお、数式１の解析解の級数和の個数は４０とした。
【００１４】
また、試料１の熱拡散率α1、試料２の熱拡散率α2を未知パラメータとし、α1を１１７±２０ｍｍ2／ｓ、α2を９６．８±２０ｍｍ2／ｓの範囲で変化させた場合のポテンシャルＪPは、図２に示すようになり多峰性が見られる。なお解析解の級数和の個数は４０とした。
【００１５】
このようにポテンシャルＪPを評価関数とした場合、極小値が見られず解が不定になる場合が予想される。また、解析解の級数和の取り方で評価関数の景観が異なり、極小値に多峰性が見られ、熱物性の同定法に滑降シンプレックス法だけでは対応できなくなる、という問題が生じる。
【００１６】
本発明は、上記の課題を解決するために創案されたものであり、積層試料を構成する物質の全ての熱物性値が未知の場合にも、熱物性値を同定することができる非定常加熱による多層積層材の熱物性同定方法を提供することを目的とする。
【課題を解決するための手段】
【００１７】
請求項１に係る発明の非定常加熱による多層積層材の熱物性同定方法は、フラッシュ加熱により得られる試料裏面の温度上昇信号と、試料の熱物性をパラメータとするフラッシュ加熱の解析解との偏差、及び、ステップ加熱により得られる試料裏面の温度上昇信号と、試料の熱物性をパラメータとするステップ加熱の解析解との偏差を評価関数とし、この評価関数が最小となるパラメータから多層積層材の各層の熱物性を決定することを特徴とする。
【００１８】
また、請求項２に係る発明の非定常加熱による多層積層材の熱物性同定方法は、請求項１に係る発明の非定常加熱による多層積層材の熱物性同定方法において、評価関数の極値探索アルゴリズムに滑降シンプレックス法を用い、評価関数の極値が多峰性を有することに対応する為、各層の熱物性の解候補である複数の粒子を含んだ粒子分布を生成する粒子フィルタを用いたことを特徴とする。なお、粒子フィルタは、予測モデルを複数用意し、観測データと尤度を比較し、適切な予測モデルを選択することで精度の高い状態を同定するものである。
【発明の効果】
【００１９】
本発明は、フラッシュ加熱の情報だけでは解が不定となることから、付帯条件としてさらにステップ加熱の情報を加えたものである。すなわち、４つの熱物性値をパラメータとした評価関数（ポテンシャル）が単一の極小値をもつなら滑降シンプレックス法で、多峰性の極値を持つなら確率要素を加味し、例えば、粒子フィルタなどの手法を取り入れ、観測された昇温曲線に解析解曲線群を評価関数が最小になるようにフィッティングさせることにより、パラメータすなわち熱物性値を決定することができる。
【図面の簡単な説明】
【００２０】
【図１】未知パラメータとして試料１の熱伝導率λ１、試料２の熱伝導率λ２を選んだ場合のポテンシャルＪＰを示すグラフである。
【図２】未知パラメータとして試料１の熱拡散率α１、試料２の熱拡散率α２を選んだ場合のポテンシャルＪＰを示すグラフである。
【図３】本発明の非定常加熱による多層積層材の熱物性同定方法を実施する装置の構成を示す図である。
【図４】フラッシュ加熱を行った場合の試料裏面の温度上昇を測定した観測信号データ及び解析解の一例を示す図である。
【図５】ステップ加熱を行った場合の試料裏面の温度上昇を測定した観測信号データ及び解析解の一例を示す図である。
【図６】未知パラメータとして試料１の熱伝導率λ１、試料２の熱伝導率λ２を選んだ場合のポテンシャルＪを示すグラフである。
【図７】未知パラメータとして試料１の熱拡散率α１、試料２の熱拡散率α２を選んだ場合のポテンシャルＪを示すグラフである。
【図８】本発明の多層積層材の熱物性同定方法を実施するアルゴリズムを示すフローチャートである。
【実施例】
【００２１】
図３は本発明の非定常加熱による多層積層材の熱物性同定方法を実施する装置の構成を示す図であり、フラッシュ加熱またはステップ加熱される被測定試料１０と、試料１０裏面の温度を測定する放射温度計２０により構成されている。図４に示すＹP(t)は、レーザパルス発生装置（図示せず）からパルス幅０．５ｍｓ程度のレーザ光で照射・加熱し、試料裏面の温度上昇を測定した場合の昇温曲線の一例であり、また、図５に示すＹS(t)は、試料１０をステップ加熱し、試料裏面の温度上昇を測定した場合の昇温曲線の一例である。
【００２２】
以下、２層積層材試料の熱物性同定を例に、フラッシュ加熱法及びステップ加熱法により得られる試料裏面の昇温信号の情報により、試料の熱伝導率、熱拡散率などの熱物性値を同定する場合について説明する。
２層積層材のフラッシュ加熱時の試料裏面温度上昇信号の解析解は、上記のように数式１で表わされるが、２層積層材のステップ加熱による裏面温度上昇信号の解析解は、下記数式３で表わされる。ただし、Ｑはステップ加熱の大きさであり、ｘkは下記数式４のｋ番目の解である。
【００２３】
【数３】
【００２４】
【数４】
【００２５】
次に、２層積層材試料を構成する各試料のそれぞれの形状（厚み）を既知として、未知の熱物性値の組み合わせについて最小２乗法で熱物性値（パラメータ）を決定する方法について説明する。
すなわち、４つの熱物性値をパラメータとした評価関数（ポテンシャル）Ｊを以下の通り定義し、このポテンシャルＪが単一の極小値を持つなら、滑降シンプレックス法で、多峰性の極値をもつなら、確率要素を加味し、例えば、粒子フィルタなどの手法を取り入れ、観測された昇温曲線に解析解曲線を評価関数が最小になるようにフィッティングさせることにより、パラメータ、すなわち、熱物性値を決定することができる。
Ｊ＝ＷSＪS＋ＷPＪP
なお、ＷS、ＷPは重み係数である。
【００２６】
上記ＪSは下記の式で表わされる。
ＪS (λ1,α1,λ2,α2)＝∫{ＹS(t)−ｙS(t, λ1,α1,λ2,α2)}２dt
なお、ＹS（ｔ）はステップ加熱した試料裏面の昇温曲線観測データ、ｙS（ｔ, λ1,α1,λ2,α2）は上記数式３で表わしたステップ加熱の場合の解析解で、図５のｙS(t)で示すものであり、ＪSはこの昇温曲線観測データＹS(t)と解析解ｙS(t)との２乗誤差の積分値である。
また、ＪＰは下記の式で表わされる。
ＪP (λ1,α1,λ2,α2)＝∫{ＹP(t) −ｙP (t, λ1,α1,λ2,α2)}２dt
なお、上記で説明したように、ＹP(t)はフラッシュ加熱した試料裏面の昇温曲線観測データ、ｙP（ｔ, λ1,α1,λ2,α2）は上記数式１で表わされるフラッシュ加熱の場合の解析解で、図４のｙP(t)で示されるものであり、ＪPはこの昇温曲線観測データＹP(t)と解析解ｙP(t)との２乗誤差の積分値である。
【００２７】
ここで、第１層目を銅（λ＝３９８Ｗ／ｍＫ、α＝１１７ｍｍ２／ｓ、厚み１ｍｍ）、第２層目をアルミ（λ＝２３７Ｗ／ｍＫ、α＝９６．８ｍｍ２／ｓ、厚み１ｍｍ）とし、銅、アルミの真値を用いて上記数式３、数式１で計算した解析解を観測データＹS (t)、ＹP(t)として用い、λ1を３９８±５０Ｗ／ｍＫ、λ2を２３７±５０Ｗ／ｍＫの範囲で変化させた場合、評価関数Ｊは、図６に示すようになり、多峰性が認められるが、極値が表れている。
【００２８】
また、試料１の熱拡散率α1と試料２の熱拡散率α2をパラメータとした評価関数Ｊは、α1を１１７±２０ｍｍ2／ｓ、α2を９６．８±２０ｍｍ2／ｓの範囲で変化させた場合、図７に示すようになり、ポテンシャルの多峰性は消えていくようである。
【００２９】
上記のように、フラッシュ加熱だけではなく、さらにステップ加熱の信号も含めると、熱伝導率の不定さが解消されることが分かり、また、極値の多峰性は探索アルゴリズムに確率的要素を持ち込めば解消できると考えられ、以下、実際の測定で熱物性を同定する場合のアルゴリズムについて、探索アルゴリズムとして滑降シンプレックス法に粒子フィルタの手法をプラスしたアルゴリズムを用いた例について、図８のフローチャートにより説明する。
【００３０】
図３の装置で、レーザパルス発生装置からパルス幅０．５ｍｓ程度のレーザ光で試料１０を照射・加熱し、放射温度計２０により試料裏面の温度上昇を測定した昇温曲線観測データＹP(t)を得、次に、試料１０をステップ加熱し、試料裏面の温度上昇を測定した昇温曲線観測データＹS (t)を得た後、図８のフローチャートのプログラムを実行する。
【００３１】
図８のフローチャートのプログラムを開始すると、まず、２層積層材の熱物性の初期値として適当な状態ベクトルｘ0を設定する（ステップ１０１）。なお、ここで言う状態ベクトルとはλ1、α1、λ2、α2を表す４元ベクトルである。
この後、昇温曲線観測データと解析解より状態ベクトルｘ0での評価関数（ポテンシャル）Ｊを計算し、Ｊ0として記憶する（ステップ１０２）。
【００３２】
次に、状態ベクトルｘ0を中心にある分散を持ったｍ個の粒子を一様に撒く（ステップ１０３）。この各粒子は色々な状態ベクトルを表している。
この後、夫々の粒子の状態ベクトル及び観測データを用いて評価関数Ｊを計算し、夫々の粒子の評価関数の極値を滑降シンプレックス法により追跡する（ステップ１０４）。なお、滑降シンプレックス法は、多次元関数の局所的な最小値を求める手法として広く用いられている方法であり、シンプレックス（Ｎ次元で説明するならば、Ｎ＋１個の頂点とそれらを結ぶ辺、面からなる多面体）を順次更新していくことで、多次元関数の局所的な最小値を求めることができる。
【００３３】
夫々の粒子が評価関数の極値に到達した後、それぞれの粒子に極値（ポテンシャル）の逆数により重み付けを行った（ステップ１０５）後、重み付けしたすべての粒子の重心を求め、新たな状態ベクトルｘとする（ステップ１０６）。
次に、新たな状態ベクトルｘでの評価関数Ｊ1を計算した（ステップ１０７）後、Ｊ0とＪ1の差が十分小さいか否かを判定し（ステップ１０８）、Ｊ0とＪ1の差が大きい場合には、Ｊ1をＪ0に置き換えた（ステップ１０９）後、ステップ１０３に戻って再び、状態ベクトルｘを中心にある分散を持ったｍ個の粒子を一様に撒く。
一方、ステップ１０８でＪ0とＪ1の差が十分小さいと判定した場合には、状態ベクトルが収束したと判断し、現在の状態ベクトルｘを解とし（ステップ１１０）、プログラムを終了する。
【００３４】
次に、試料１、２の厚みｄ1、ｄ2が既知で、λ1、α1、λ2、α2が未知であり、目標値で演算した解析解を昇温曲線観測データとして用い、各パラメータによる解析解から評価関数を計算して上記探索アルゴリズムにより最適解を探索した結果について説明する。
下記表１の条件で以下のソフトアルゴリズムを実行した。
１．状態ベクトルｘの初期値を与える。
２．分散SigmaＷにより状態ベクトルｘを中心に分布する粒子１００個を撒く。
３．それぞれの粒子の評価関数を計算し、滑降シンプレックス法により評価関数の極値を追跡する。
４．それぞれの粒子に到達極値の逆数により重み付けを行う。
５．全ての粒子の重心を求め、新たな状態ベクトルｘとする。
６．以下収束するまで２〜５を繰り返す。
７．収束判定後、最小値を与える状態ベクトルを解とする。
【００３５】
【表１】
【００３６】
上記探索アルゴリズムによる最適解の探索結果は下記表２に示す通りであり、初期値にもよると思われるが、３回の試行で正解にたどり着いている。なお、表２中の最小値Ｘは、計算実行中に現れた極値の最小値を示した状態ベクトルの解であり、このように最小値を示した状態ベクトルを求めておくことにより、初期設定の状態ベクトルから離れた位置に解がある場合にも対応することができる。
【００３７】
【表２】
【００３８】
また、実際によくある例として、試料１の熱物性値及び試料全体の厚みを既知として、試料２の厚みｄ2、熱伝導率λ2、熱拡散率α2が未知の場合に、下記表３の条件で以下のソフトアルゴリズムを実行した。
１．状態ベクトルｘの初期値を与える。
２．分散SigmaＷにより状態ベクトルｘを中心に分布する粒子１０個を撒く。
３．それぞれの粒子の評価関数を計算し、滑降シンプレックス法により評価関数の極値を追跡する。
４．それぞれの粒子に到達極値の逆数により重み付けを行う。
５．全ての粒子の重心を求め、新たな状態ベクトルｘとする。
６．以下収束するまで２〜５を繰り返す。
７．収束判定後、最小値を与える状態ベクトルを解とする。
【００３９】
【表３】
【００４０】
上記探索アルゴリズムによる最適解の探索結果は下記表４に示す通りであり、非常に収束性がよく停留点がなく、粒子数１０個程度で十分正解にたどり着くことができる。なお、上記と同様に、表４中の最小値Ｘは、計算実行中に現れた極値の最小値を示した状態ベクトルの解である。
【００４１】
【表４】
【００４２】
さらに、２つの試料の４つの熱物性値が不明で、試料形状は全体の厚さＬしかわからない場合を例として、下記表５の条件で以下のソフトアルゴリズムを実行した。
１．状態ベクトルｘの初期値を与える。
２．分散SigmaＷにより状態ベクトルｘを中心に分布する粒子１００個を撒く。
３．それぞれの粒子の評価関数を計算し、滑降シンプレックス法により評価関数の極値を追跡する。
４．それぞれの粒子に到達極値の逆数により重み付けを行う。
５．全ての粒子の重心を求め、新たな状態ベクトルｘとする。
６．以下収束するまで２〜５を繰り返す。
７．収束判定後、最小値を与える状態ベクトルを解とする。
【００４３】
【表５】
【００４４】
上記探索アルゴリズムによる最適解の探索結果は下記表６に示す通りであり、精度は落ちるが同定できることが分かる。なお、上記と同様に、表６中の最小値Ｘは、計算実行中に現れた極値の最小値を示した状態ベクトルの解である。
【００４５】
【表６】
【００４６】
積層材の熱物性の同定が困難であった理由の一つはポテンシャルが極値を持たないため不定になる可能性があること、もう一つは極値が単峰性でなく多峰性となることと思われるが、本発明による非定常加熱による多層積層材の熱物性同定方法によれば、不定になることについては目的関数に付加条件としてステップ加熱情報を追加することで、多峰性については最小値探索アルゴリズムに確率性を持ち込むことで解決することができた。
【００４７】
なお、上記の実施例のアルゴリズムでは、重み付けした粒子群の重心を求め、これを新たな状態ベクトルとしたが、評価関数の最小値を与える状態ベクトルを新たな状態ベクトルとすることもできる。
また、上記の実施例では、極値を探索するのに滑降シンプレックス法を採用したが、山登り法等を使用することも可能である。また、確率的要素として粒子フィルタのアルゴリズムを用いたが、これを遺伝的アルゴリズムに置き換えることも可能である。
さらに、上記の実施例では、観測データと解析解の偏差を求めるのに、２乗誤差を使用したが、偏差として誤差の絶対値の和等を使用することも可能である。
【００４８】
また、上記の実施例では、観測データとの比較に解析解を用いたが、解析解が得られていない場合は、この種の問題の解法の多くがラプラス変換をして問題を解いており、数値解析によるラプラス逆変換より得た解を利用することも可能である。
【００４９】
さらに、上記実施例では、付加条件としてステップ加熱信号を使用したが、図５に示すグラフの直線の傾きの逆数が試料の全熱容量を示すから、ステップ加熱信号を別途測定した試料の全熱容量と置き換えてもよい。但しこのときは、付加条件の情報量が少なくなるから、例えば表３、表４に示す例の様に、当然決定できる未知変数も少なくなる。
【符号の説明】
【００５０】
１０ 試料
２０ 放射温度計
【技術分野】
【０００１】
本発明は、試料の熱伝導率、熱拡散率などの熱物性の同定方法、特に、多層積層材の各層の熱伝導率、熱拡散率などの熱物性を同定する方法に関する。
【背景技術】
【０００２】
金属、セラミックス、ガラス、カーボン、プラスチックなど個体材料の熱的特性を精度よく測定する方法としてフラッシュ加熱法がある。このフラッシュ加熱法は、試料の表面にレーザ光を瞬間的（フラッシュ状）に照射して加熱し、一定距離離れた場所に伝わる熱を温度変化から測定するものである。例えば、直径１０ｍｍ、厚さ１〜３ｍｍ程度の小円板状試料の片面をレーザ光で瞬時に加熱して反対側の温度挙動から熱測定を行う装置があり、熱的特性としては、熱拡散率、及び比熱容量を測定することができ、さらに密度が既知であれば、熱伝導率を算出することができる。
【０００３】
一方、近年の新しい材料開発の動向の一つとして、複数の素材を組み合わせることにより、今までに無かった新しい機能を持たせた材料を開発することが行われており、その中でも特に高熱伝導性基板などが注目されている。しかし、この分野で必要とされる評価測定の部門で、出来上がった新機能材料の各々の層の熱物性を直接測定する方法はまだ確立されていない。
【０００４】
このような多層試料の熱物性を測定する場合、従来、定常法での熱伝導率測定により、多層試料全体の熱伝導率を測定することが行われている。この測定法は試料板の片面を加熱、反対面を冷却し、試料板厚み方向に温度勾配をつけ、その温度差と試料厚み及び試料を貫通する熱流束から熱伝導率を算出する方法である。
【０００５】
また、各層の熱物性を測定する試みとして、上記のフラッシュ加熱法を使用することが提案されている。すなわち、未知の一層を含む多層材料の表面側からレーザフラッシュを照射して得られる裏面側の温度応答と、多層材料の熱物性による理論温度応答との２乗誤差を評価関数（ポテンシャル）とし、この評価関数が最小値になるパラメータを求めることにより高精度でかつ効率的に多層材料中の一層における未知の熱拡散率、比熱等の熱定数を決定するようにしている（例えば、特許文献１参照）。
【先行技術文献】
【特許文献】
【０００６】
【特許文献１】特開平８−２６１９６７号公報
【発明の概要】
【発明が解決しようとする課題】
【０００７】
上記のように、熱拡散率及び比熱が未知の一層を含む多層材料の表面側からレーザフラッシュを照射して得られる裏面側の温度応答と、熱拡散率、比熱、及びビオー数を変数として含む理論温度応答とを比較することにより、多層積層材の熱物性を同定することが提案されているが、この方法では、多層材料の一層の熱拡散率及び比熱が未知の場合のみ、熱物性を同定できるだけであり、積層試料を構成する物質の全ての熱物性値が未知の場合には同定することは困難である。
【０００８】
以下、２層積層試料を構成する各試料をそれぞれの形状（厚み）を既知として、未知の熱物性値の組み合わせについてフラッシュ加熱法を用いて最小２乗法で各パラメータ、すなわち熱物性値を決定できるかについて説明する。
未知パラメータを独立変数とした評価関数（ポテンシャル）を以下のＪPで定義し、このポテンシャルが単一の極値を持つなら、滑降シンプレックス法で観測された昇温曲線に解析解曲線群をフィッティングさせることによりパラメータ、すなわち、熱物性値を決定することができる。
ＪP(λ1,α1,λ2,α2)=∫{ＹP(t)-ｙP(t,λ1,α1,λ2,α2)}2dt
ただし、ＹP(t)はフラッシュ加熱で得られた観測信号データ、ｙP(t,λ1,α1,λ2,α2)はフィッティングさせる解析解を表す。
【０００９】
上記の解析解ｙP (t,λ1,α1,λ2,α2)は、以下の数式１のθr（ｔ）で表わされる。なお、数式１において、λ、α、ｄをそれぞれ熱伝導率、熱拡散率、試料厚みとしたとき、ｂは熱浸透率でｂ＝λ／√α、τは特性時間でτ＝ｄ2／αである。また、添え字１、２は１層目、２層目の物性値を表す。
【００１０】
【数１】
【００１１】
ただし、ｓkはｓを複素数として下記数式２のｋ番目の解である。
【００１２】
【数２】
【００１３】
ここで、第１層目（試料１）を銅（λ＝３９８Ｗ／ｍＫ、α＝１１７ｍｍ2／ｓ、厚み１ｍｍ）、第２層目（試料２）をアルミ（λ＝２３７Ｗ／ｍＫ、α＝９６．８ｍｍ2／ｓ、厚み１ｍｍ）とし、また、昇温波形の最大値は、試料からの熱放散がないと考え（すなわちビオー数をゼロとした）、観測データＹP(t)から知ることができるとし、λ1＝３９８ｍＫ、α1＝１１７ｍｍ2／ｓ、λ2＝２３７Ｗ／ｍＫ、α2＝９６．８ｍｍ2／ｓのときの上記数式１で計算した解析解を観測データＹP(t)として用い、試料１の熱伝導率λ1、試料２の熱伝導率λ2を未知パラメータとし、λ1を３９８±５０Ｗ／ｍＫ、λ2を２３７±５０Ｗ／ｍＫの範囲で変化させた場合のポテンシャルＪPは、図１に示すようになり、極値の識別は困難である。なお、数式１の解析解の級数和の個数は４０とした。
【００１４】
また、試料１の熱拡散率α1、試料２の熱拡散率α2を未知パラメータとし、α1を１１７±２０ｍｍ2／ｓ、α2を９６．８±２０ｍｍ2／ｓの範囲で変化させた場合のポテンシャルＪPは、図２に示すようになり多峰性が見られる。なお解析解の級数和の個数は４０とした。
【００１５】
このようにポテンシャルＪPを評価関数とした場合、極小値が見られず解が不定になる場合が予想される。また、解析解の級数和の取り方で評価関数の景観が異なり、極小値に多峰性が見られ、熱物性の同定法に滑降シンプレックス法だけでは対応できなくなる、という問題が生じる。
【００１６】
本発明は、上記の課題を解決するために創案されたものであり、積層試料を構成する物質の全ての熱物性値が未知の場合にも、熱物性値を同定することができる非定常加熱による多層積層材の熱物性同定方法を提供することを目的とする。
【課題を解決するための手段】
【００１７】
請求項１に係る発明の非定常加熱による多層積層材の熱物性同定方法は、フラッシュ加熱により得られる試料裏面の温度上昇信号と、試料の熱物性をパラメータとするフラッシュ加熱の解析解との偏差、及び、ステップ加熱により得られる試料裏面の温度上昇信号と、試料の熱物性をパラメータとするステップ加熱の解析解との偏差を評価関数とし、この評価関数が最小となるパラメータから多層積層材の各層の熱物性を決定することを特徴とする。
【００１８】
また、請求項２に係る発明の非定常加熱による多層積層材の熱物性同定方法は、請求項１に係る発明の非定常加熱による多層積層材の熱物性同定方法において、評価関数の極値探索アルゴリズムに滑降シンプレックス法を用い、評価関数の極値が多峰性を有することに対応する為、各層の熱物性の解候補である複数の粒子を含んだ粒子分布を生成する粒子フィルタを用いたことを特徴とする。なお、粒子フィルタは、予測モデルを複数用意し、観測データと尤度を比較し、適切な予測モデルを選択することで精度の高い状態を同定するものである。
【発明の効果】
【００１９】
本発明は、フラッシュ加熱の情報だけでは解が不定となることから、付帯条件としてさらにステップ加熱の情報を加えたものである。すなわち、４つの熱物性値をパラメータとした評価関数（ポテンシャル）が単一の極小値をもつなら滑降シンプレックス法で、多峰性の極値を持つなら確率要素を加味し、例えば、粒子フィルタなどの手法を取り入れ、観測された昇温曲線に解析解曲線群を評価関数が最小になるようにフィッティングさせることにより、パラメータすなわち熱物性値を決定することができる。
【図面の簡単な説明】
【００２０】
【図１】未知パラメータとして試料１の熱伝導率λ１、試料２の熱伝導率λ２を選んだ場合のポテンシャルＪＰを示すグラフである。
【図２】未知パラメータとして試料１の熱拡散率α１、試料２の熱拡散率α２を選んだ場合のポテンシャルＪＰを示すグラフである。
【図３】本発明の非定常加熱による多層積層材の熱物性同定方法を実施する装置の構成を示す図である。
【図４】フラッシュ加熱を行った場合の試料裏面の温度上昇を測定した観測信号データ及び解析解の一例を示す図である。
【図５】ステップ加熱を行った場合の試料裏面の温度上昇を測定した観測信号データ及び解析解の一例を示す図である。
【図６】未知パラメータとして試料１の熱伝導率λ１、試料２の熱伝導率λ２を選んだ場合のポテンシャルＪを示すグラフである。
【図７】未知パラメータとして試料１の熱拡散率α１、試料２の熱拡散率α２を選んだ場合のポテンシャルＪを示すグラフである。
【図８】本発明の多層積層材の熱物性同定方法を実施するアルゴリズムを示すフローチャートである。
【実施例】
【００２１】
図３は本発明の非定常加熱による多層積層材の熱物性同定方法を実施する装置の構成を示す図であり、フラッシュ加熱またはステップ加熱される被測定試料１０と、試料１０裏面の温度を測定する放射温度計２０により構成されている。図４に示すＹP(t)は、レーザパルス発生装置（図示せず）からパルス幅０．５ｍｓ程度のレーザ光で照射・加熱し、試料裏面の温度上昇を測定した場合の昇温曲線の一例であり、また、図５に示すＹS(t)は、試料１０をステップ加熱し、試料裏面の温度上昇を測定した場合の昇温曲線の一例である。
【００２２】
以下、２層積層材試料の熱物性同定を例に、フラッシュ加熱法及びステップ加熱法により得られる試料裏面の昇温信号の情報により、試料の熱伝導率、熱拡散率などの熱物性値を同定する場合について説明する。
２層積層材のフラッシュ加熱時の試料裏面温度上昇信号の解析解は、上記のように数式１で表わされるが、２層積層材のステップ加熱による裏面温度上昇信号の解析解は、下記数式３で表わされる。ただし、Ｑはステップ加熱の大きさであり、ｘkは下記数式４のｋ番目の解である。
【００２３】
【数３】
【００２４】
【数４】
【００２５】
次に、２層積層材試料を構成する各試料のそれぞれの形状（厚み）を既知として、未知の熱物性値の組み合わせについて最小２乗法で熱物性値（パラメータ）を決定する方法について説明する。
すなわち、４つの熱物性値をパラメータとした評価関数（ポテンシャル）Ｊを以下の通り定義し、このポテンシャルＪが単一の極小値を持つなら、滑降シンプレックス法で、多峰性の極値をもつなら、確率要素を加味し、例えば、粒子フィルタなどの手法を取り入れ、観測された昇温曲線に解析解曲線を評価関数が最小になるようにフィッティングさせることにより、パラメータ、すなわち、熱物性値を決定することができる。
Ｊ＝ＷSＪS＋ＷPＪP
なお、ＷS、ＷPは重み係数である。
【００２６】
上記ＪSは下記の式で表わされる。
ＪS (λ1,α1,λ2,α2)＝∫{ＹS(t)−ｙS(t, λ1,α1,λ2,α2)}２dt
なお、ＹS（ｔ）はステップ加熱した試料裏面の昇温曲線観測データ、ｙS（ｔ, λ1,α1,λ2,α2）は上記数式３で表わしたステップ加熱の場合の解析解で、図５のｙS(t)で示すものであり、ＪSはこの昇温曲線観測データＹS(t)と解析解ｙS(t)との２乗誤差の積分値である。
また、ＪＰは下記の式で表わされる。
ＪP (λ1,α1,λ2,α2)＝∫{ＹP(t) −ｙP (t, λ1,α1,λ2,α2)}２dt
なお、上記で説明したように、ＹP(t)はフラッシュ加熱した試料裏面の昇温曲線観測データ、ｙP（ｔ, λ1,α1,λ2,α2）は上記数式１で表わされるフラッシュ加熱の場合の解析解で、図４のｙP(t)で示されるものであり、ＪPはこの昇温曲線観測データＹP(t)と解析解ｙP(t)との２乗誤差の積分値である。
【００２７】
ここで、第１層目を銅（λ＝３９８Ｗ／ｍＫ、α＝１１７ｍｍ２／ｓ、厚み１ｍｍ）、第２層目をアルミ（λ＝２３７Ｗ／ｍＫ、α＝９６．８ｍｍ２／ｓ、厚み１ｍｍ）とし、銅、アルミの真値を用いて上記数式３、数式１で計算した解析解を観測データＹS (t)、ＹP(t)として用い、λ1を３９８±５０Ｗ／ｍＫ、λ2を２３７±５０Ｗ／ｍＫの範囲で変化させた場合、評価関数Ｊは、図６に示すようになり、多峰性が認められるが、極値が表れている。
【００２８】
また、試料１の熱拡散率α1と試料２の熱拡散率α2をパラメータとした評価関数Ｊは、α1を１１７±２０ｍｍ2／ｓ、α2を９６．８±２０ｍｍ2／ｓの範囲で変化させた場合、図７に示すようになり、ポテンシャルの多峰性は消えていくようである。
【００２９】
上記のように、フラッシュ加熱だけではなく、さらにステップ加熱の信号も含めると、熱伝導率の不定さが解消されることが分かり、また、極値の多峰性は探索アルゴリズムに確率的要素を持ち込めば解消できると考えられ、以下、実際の測定で熱物性を同定する場合のアルゴリズムについて、探索アルゴリズムとして滑降シンプレックス法に粒子フィルタの手法をプラスしたアルゴリズムを用いた例について、図８のフローチャートにより説明する。
【００３０】
図３の装置で、レーザパルス発生装置からパルス幅０．５ｍｓ程度のレーザ光で試料１０を照射・加熱し、放射温度計２０により試料裏面の温度上昇を測定した昇温曲線観測データＹP(t)を得、次に、試料１０をステップ加熱し、試料裏面の温度上昇を測定した昇温曲線観測データＹS (t)を得た後、図８のフローチャートのプログラムを実行する。
【００３１】
図８のフローチャートのプログラムを開始すると、まず、２層積層材の熱物性の初期値として適当な状態ベクトルｘ0を設定する（ステップ１０１）。なお、ここで言う状態ベクトルとはλ1、α1、λ2、α2を表す４元ベクトルである。
この後、昇温曲線観測データと解析解より状態ベクトルｘ0での評価関数（ポテンシャル）Ｊを計算し、Ｊ0として記憶する（ステップ１０２）。
【００３２】
次に、状態ベクトルｘ0を中心にある分散を持ったｍ個の粒子を一様に撒く（ステップ１０３）。この各粒子は色々な状態ベクトルを表している。
この後、夫々の粒子の状態ベクトル及び観測データを用いて評価関数Ｊを計算し、夫々の粒子の評価関数の極値を滑降シンプレックス法により追跡する（ステップ１０４）。なお、滑降シンプレックス法は、多次元関数の局所的な最小値を求める手法として広く用いられている方法であり、シンプレックス（Ｎ次元で説明するならば、Ｎ＋１個の頂点とそれらを結ぶ辺、面からなる多面体）を順次更新していくことで、多次元関数の局所的な最小値を求めることができる。
【００３３】
夫々の粒子が評価関数の極値に到達した後、それぞれの粒子に極値（ポテンシャル）の逆数により重み付けを行った（ステップ１０５）後、重み付けしたすべての粒子の重心を求め、新たな状態ベクトルｘとする（ステップ１０６）。
次に、新たな状態ベクトルｘでの評価関数Ｊ1を計算した（ステップ１０７）後、Ｊ0とＪ1の差が十分小さいか否かを判定し（ステップ１０８）、Ｊ0とＪ1の差が大きい場合には、Ｊ1をＪ0に置き換えた（ステップ１０９）後、ステップ１０３に戻って再び、状態ベクトルｘを中心にある分散を持ったｍ個の粒子を一様に撒く。
一方、ステップ１０８でＪ0とＪ1の差が十分小さいと判定した場合には、状態ベクトルが収束したと判断し、現在の状態ベクトルｘを解とし（ステップ１１０）、プログラムを終了する。
【００３４】
次に、試料１、２の厚みｄ1、ｄ2が既知で、λ1、α1、λ2、α2が未知であり、目標値で演算した解析解を昇温曲線観測データとして用い、各パラメータによる解析解から評価関数を計算して上記探索アルゴリズムにより最適解を探索した結果について説明する。
下記表１の条件で以下のソフトアルゴリズムを実行した。
１．状態ベクトルｘの初期値を与える。
２．分散SigmaＷにより状態ベクトルｘを中心に分布する粒子１００個を撒く。
３．それぞれの粒子の評価関数を計算し、滑降シンプレックス法により評価関数の極値を追跡する。
４．それぞれの粒子に到達極値の逆数により重み付けを行う。
５．全ての粒子の重心を求め、新たな状態ベクトルｘとする。
６．以下収束するまで２〜５を繰り返す。
７．収束判定後、最小値を与える状態ベクトルを解とする。
【００３５】
【表１】
【００３６】
上記探索アルゴリズムによる最適解の探索結果は下記表２に示す通りであり、初期値にもよると思われるが、３回の試行で正解にたどり着いている。なお、表２中の最小値Ｘは、計算実行中に現れた極値の最小値を示した状態ベクトルの解であり、このように最小値を示した状態ベクトルを求めておくことにより、初期設定の状態ベクトルから離れた位置に解がある場合にも対応することができる。
【００３７】
【表２】
【００３８】
また、実際によくある例として、試料１の熱物性値及び試料全体の厚みを既知として、試料２の厚みｄ2、熱伝導率λ2、熱拡散率α2が未知の場合に、下記表３の条件で以下のソフトアルゴリズムを実行した。
１．状態ベクトルｘの初期値を与える。
２．分散SigmaＷにより状態ベクトルｘを中心に分布する粒子１０個を撒く。
３．それぞれの粒子の評価関数を計算し、滑降シンプレックス法により評価関数の極値を追跡する。
４．それぞれの粒子に到達極値の逆数により重み付けを行う。
５．全ての粒子の重心を求め、新たな状態ベクトルｘとする。
６．以下収束するまで２〜５を繰り返す。
７．収束判定後、最小値を与える状態ベクトルを解とする。
【００３９】
【表３】
【００４０】
上記探索アルゴリズムによる最適解の探索結果は下記表４に示す通りであり、非常に収束性がよく停留点がなく、粒子数１０個程度で十分正解にたどり着くことができる。なお、上記と同様に、表４中の最小値Ｘは、計算実行中に現れた極値の最小値を示した状態ベクトルの解である。
【００４１】
【表４】
【００４２】
さらに、２つの試料の４つの熱物性値が不明で、試料形状は全体の厚さＬしかわからない場合を例として、下記表５の条件で以下のソフトアルゴリズムを実行した。
１．状態ベクトルｘの初期値を与える。
２．分散SigmaＷにより状態ベクトルｘを中心に分布する粒子１００個を撒く。
３．それぞれの粒子の評価関数を計算し、滑降シンプレックス法により評価関数の極値を追跡する。
４．それぞれの粒子に到達極値の逆数により重み付けを行う。
５．全ての粒子の重心を求め、新たな状態ベクトルｘとする。
６．以下収束するまで２〜５を繰り返す。
７．収束判定後、最小値を与える状態ベクトルを解とする。
【００４３】
【表５】
【００４４】
上記探索アルゴリズムによる最適解の探索結果は下記表６に示す通りであり、精度は落ちるが同定できることが分かる。なお、上記と同様に、表６中の最小値Ｘは、計算実行中に現れた極値の最小値を示した状態ベクトルの解である。
【００４５】
【表６】
【００４６】
積層材の熱物性の同定が困難であった理由の一つはポテンシャルが極値を持たないため不定になる可能性があること、もう一つは極値が単峰性でなく多峰性となることと思われるが、本発明による非定常加熱による多層積層材の熱物性同定方法によれば、不定になることについては目的関数に付加条件としてステップ加熱情報を追加することで、多峰性については最小値探索アルゴリズムに確率性を持ち込むことで解決することができた。
【００４７】
なお、上記の実施例のアルゴリズムでは、重み付けした粒子群の重心を求め、これを新たな状態ベクトルとしたが、評価関数の最小値を与える状態ベクトルを新たな状態ベクトルとすることもできる。
また、上記の実施例では、極値を探索するのに滑降シンプレックス法を採用したが、山登り法等を使用することも可能である。また、確率的要素として粒子フィルタのアルゴリズムを用いたが、これを遺伝的アルゴリズムに置き換えることも可能である。
さらに、上記の実施例では、観測データと解析解の偏差を求めるのに、２乗誤差を使用したが、偏差として誤差の絶対値の和等を使用することも可能である。
【００４８】
また、上記の実施例では、観測データとの比較に解析解を用いたが、解析解が得られていない場合は、この種の問題の解法の多くがラプラス変換をして問題を解いており、数値解析によるラプラス逆変換より得た解を利用することも可能である。
【００４９】
さらに、上記実施例では、付加条件としてステップ加熱信号を使用したが、図５に示すグラフの直線の傾きの逆数が試料の全熱容量を示すから、ステップ加熱信号を別途測定した試料の全熱容量と置き換えてもよい。但しこのときは、付加条件の情報量が少なくなるから、例えば表３、表４に示す例の様に、当然決定できる未知変数も少なくなる。
【符号の説明】
【００５０】
１０ 試料
２０ 放射温度計
【特許請求の範囲】
【請求項１】
フラッシュ加熱により得られる試料裏面の温度上昇信号と、試料の熱物性をパラメータとするフラッシュ加熱の解析解との偏差、及び、ステップ加熱により得られる試料裏面の温度上昇信号と、試料の熱物性をパラメータとするステップ加熱の解析解との偏差を評価関数とし、この評価関数が最小となるパラメータから多層積層材の各層の熱物性を決定することを特徴とする、非定常加熱による多層積層材の熱物性同定方法。
【請求項２】
請求項１に係る発明の非定常加熱による多層積層材の熱物性同定方法において、評価関数の極値探索アルゴリズムに滑降シンプレックス法を用い、評価関数の極値が多峰性を有することに対応する為、各層の熱物性の解候補である複数の粒子を含んだ粒子分布を生成する粒子フィルタを用いたことを特徴とする、非定常加熱による多層積層材の熱物性同定方法。
【請求項１】
フラッシュ加熱により得られる試料裏面の温度上昇信号と、試料の熱物性をパラメータとするフラッシュ加熱の解析解との偏差、及び、ステップ加熱により得られる試料裏面の温度上昇信号と、試料の熱物性をパラメータとするステップ加熱の解析解との偏差を評価関数とし、この評価関数が最小となるパラメータから多層積層材の各層の熱物性を決定することを特徴とする、非定常加熱による多層積層材の熱物性同定方法。
【請求項２】
請求項１に係る発明の非定常加熱による多層積層材の熱物性同定方法において、評価関数の極値探索アルゴリズムに滑降シンプレックス法を用い、評価関数の極値が多峰性を有することに対応する為、各層の熱物性の解候補である複数の粒子を含んだ粒子分布を生成する粒子フィルタを用いたことを特徴とする、非定常加熱による多層積層材の熱物性同定方法。
【図１】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【公開番号】特開２０１２−６８１９３（Ｐ２０１２−６８１９３Ａ）
【公開日】平成２４年４月５日（２０１２．４．５）
【国際特許分類】
【出願番号】特願２０１０−２１５０１５（Ｐ２０１０−２１５０１５）
【出願日】平成２２年９月２７日（２０１０．９．２７）
【出願人】（０００１６１９３２）京都電子工業株式会社 (29)
【Ｆターム（参考）】
【公開日】平成２４年４月５日（２０１２．４．５）
【国際特許分類】
【出願日】平成２２年９月２７日（２０１０．９．２７）
【出願人】（０００１６１９３２）京都電子工業株式会社 (29)
【Ｆターム（参考）】
[ Back to top ]