長岡工業高等専門学校研究紀要

論第 52 巻（2016）

文

非破壊打音試験における打音信号同定のための逆解析
倉橋 貴彦1・松岡 瑛喜2・村上 祐貴3・外山 茂浩4・池田 富士雄5・井山 徹郎6・田原 喜宏7
1 長岡技術科学大学大学院（Department of Mechanical Engineering, Nagaoka University of Technology）
2 長岡技術科学大学大学院（Department of Mechanical Engineering, Nagaoka University of Technology）
3 環境都市工学科（Department of Civil Engineering, National Institute of Technology, Nagaoka College）
4 電子制御工学科（Department of Electronic Control Engineering,

National Institute of Technology, Nagaoka College）
5 機械工学科（Department ofMechanical Engineering, National Institute of Technology, Nagaoka College）

6 機械工学科（Department ofMechanical Engineering, National Institute of Technology, Nagaoka College）

7 一般教育科－数学（Liberal Arts-Mathematics, National Institute of Technology, Nagaoka College）

Inverse analysis for Identification of Hammering Signal
in Non-destructive Hammering Test
Takahiko KURAHASHI1, Eiki MATSUOKA2, Yuki MURAKAMI3,
Shigehiro TOYAMA4, Fujio IKEDA5,Tetsuro IYAMA6 and Yoshihiro TAWARA7Abstract
In this study, we present the inverse analysis for identification of hammering signal in non-destructive
hammering test. As the observation value, the observed voltage from microphone is transformed to soundpressure level, and the performance function is defined by square sum of residual between the obtained
and the computed sound pressure levels. Here, the problem is to find the input sound pressure level so as
to minimize the performance function. The formulation for this problem is carried out by the adjoint
variable method, and the numerical simulation of the sound pressure level propagation is carried out
based on the wave equation and the finite element method.
Key Words : Non-destructive testing, Hammering test, Wave equation, Finite element method, Adjoint
variable method

１．序論

おいても，構造物の維持管理の検査を行う際，様々な箇所で適用されている．経験則に基づき，検査者

非破壊試験法の一つである打音試験法は，現在にの判断で構造物の不良箇所を特定することもあるが，- 34 -

倉橋 貴彦・松岡 瑛喜・村上 祐貴・外山 茂浩・池田 富士雄・井山 徹郎・田原 喜宏

マイクロフォン等で，打音時の音を計測し，構造物内の欠陥を探査する研究も行われている 1)．論文内

p  unknown on 1(5)b  c 2 p ,i ni  bˆ on  2(6)では，打撃時における変位波形と，マイクロフォンによって測定された電圧の比較をしており，実験と解析の値は一致していないものの，モデルの改善により解析による欠陥位置等の事前推定が可能であることも言及している．そこで，本論文では，打撃にここで c は音速，ni は外向きの単位法線ベクトル，よる変位波形でなく，電圧波形から音圧レベルを算

ハットは既知量，Γ1 は音圧を入力する境界（打音出し，音圧の伝播計算をすることで，同等の物理量点），Γ2 はそれ以外の境界を意味する．

を比較することで，実験と解析の結果を行うことを考える．この検討の場合，解析における境界条件の2.2 停留条件

設定として，打撃時における音圧波形を知る必要が式(2)～式(6)に示した拘束条件のもと評価関数をある．そこで，本検討では，マイクロフォンにより最小とする問題を解く際，随伴変数 λ を導入するこ測定した電圧信号を音圧レベルに換算し，得られたとにより式(7)に示すラグランジュ関数を得る．

データから打撃点における音圧レベルを同定することを試みる．音圧レベルの波形の伝播の解析には有J  J 

限要素法を適用し，打撃点における音圧レベルの同   p  c p ddt
tft02,ii

(7)

定には，随伴変数法を適用する．
ここで，ラグランジュ関数に対する停留条件の導出を行う．ラグランジュ関数 J*の第一変分を計算す

２．入力信号同定のための計算アルゴリズム

ることにより随伴方程式および随伴変数に対する終端・境界条件として式(8)～式(12)が導かれる．

2.1 入力信号同定のための評価関数の定義
本論文では，マイクロフォンにより計測された電圧値を音圧レベルに換算し，基準音圧値を音圧レベルの式に代入し得られる音圧値を評価関数に用いる．
計測電圧から求められた音圧を pobs.とし，シミュレ  c 2 ,ii   p  pobs. Q  0(8)  0 at t  t fin (9)  0 at t  t fin (10)ーションにより求まる音圧値 p の残差二乗和を評価関数（式(1)）として定義する．

J  p  p

1 tf
2 t0







obs. Q p  p obs. ddt

  0 on 1(1)s  c 2 ,i ni  0

ここに tf はシミュレーションの終端時刻， Q は重on 2(11)(12)

み定数であり，重み定数 Q は計測点において 1，それ以外においては 0 と設定する．本論文における入2.3 状態方程式と随伴方程式の離散化

力電圧同定問題においては，式(2)に示す波動方程

状態変数および随伴変数に対して，それぞれ速度式および式(3)～式(6)に示す初期・境界条件のもとの変数 v と μ を導入し，重み関数を乗じ，各要素領で評価関数を最小とする打音点における音圧の経時域において積分することにより，状態変数に対する変化を求めることになる．

重み付き残差方程式（式(13)，(14)）および随伴変数に対する重み付き残差方程式（式(15)，(16)）がp  c 2 p,ii  0(2)p  pˆ at t  t 0in (3)p  p̂ at t  t 0in (4)得られる．

- 35 - v (v  c p )d  0e2,ii(13)非破壊打音試験における打音信号同定のための逆解析e

p ( p  v)d  0

    c   Q( p  p )d  0
2,iieobs.

If J(l+1)≦J(l) then α(l+1)=0.5α(l)goto Step3(14)(15)

Else α(l+1)=1.1α(l) go to Step4

３．数値実験
3.1 数値解析例１

  (   )d  0e

数値実験の検討モデル図を図-1に示す．音圧の入(16)力点を(x,y)=(0.5m,1.0m)の箇所に設定し，音圧の計

ここで，空間方向の離散化に対しては三角形一次の形状関数を用いたガラーキン法を適用し，また時間方向の離散化に対しては Krenk の提案する手法 3)を適用することで最終的に式(13)～式(16)に対する有限要素方程式を得ることになる．以上の定式化により導出された有限要素方程式の計算を行うことにより，計算領域内における状態変数および随伴変数の値を数値的に得ることができる．

測点を(x,y)=(0.8m,1.0m)に設定する．音圧を加える点以外の外部境界においては，音圧に関する流束ベクトル(式(6))を零と設定する．ここで検討内容としては，計測点における音圧を計測音圧に近づけるような，入力点における音圧波形を求める問題について逆解析を行う．打音試験における入力点での音圧は，パルス入力により表されるが，入力点における波形が本論文に示す方法により経時的に正しく求まるのか確認を行うため，本検討では，式(18)に示す正弦波により表される音圧波形を，入力点におけ2.4 ラグランジュ関数の最小化手法ラグランジュ関数の最小化をするために，最急降る正解値として仮定する．

下法 4)を導入する．打音点における音圧の更新式は
 2 p  sint T 式(17)のように書くことができる．

p (l 1)  p (l )   (l )J (l )p (l )

on 1(18)

式(18)において，Tは周期を示し，T=100sec.と設(17)定する．正解となる音圧（式(18)により与えられた音圧）を入力点に与えた順解析を行い，計測点におここに l は反復回数を示し，α はステップ長さを示いて得られた波形を人工的な計測値として使用する．

す．ここで最急降下法に基づく計算アルゴリズムは本検討では，入力点において，音圧が経時的に零と以下のようになる．

いう条件から逆解析の計算を始め，前章に示した計算アルゴリズムに従い反復計算を行い，計測値からStep1 境界 Γ1 上における打音点における音圧の初期値 p(l)および収束判定定数 ε を設定する．
Step2 状態方程式より p(l)を算定し，評価関数 J(l)を式(18)に示す音圧波形が入力点において求まるか検討を行う．ここに，計算条件を表-1に示す．
Input point of sound pressure

計算する．

0.8m

Step3 随伴方程式より λ(l)を算定し，ラグランジュ関

0.5m

Observation point

数 J*(l)の音圧 p(l)に対する勾配∂J*(l)/∂ p(l)を計算
する．
Step4 式(18)をもとに音圧を更新し p(l+1)とする．
Step5 状態方程式より p(l+1)を算定し，評価関数 J(l+1)
を計算する．
Step6 収束判定を行う．
If |J(l+1)－J(l)|<ε then Stop
Else go to Step 7

図-1

Step7 ステップ長さ α(l)の更新を行う．

- 36 -計算モデル図

倉橋 貴彦・松岡 瑛喜・村上 祐貴・外山 茂浩・池田 富士雄・井山 徹郎・田原 喜宏

表-1 計算条件

総節点数
総要素数
時間増分量Δt，sec.
波速，m/ sec.
時間ステップ
Sakawa-Shindo法における1イタレーション目の重み関数の対角値
1イタレーション目の入力点における音圧，Pa.
収束判定定数ε2601
5000
0.04
0.01
10000
1.00.010-6
図-2

評価関数の収束履歴

3.2 入力信号同定の逆解析
前節に示す計算条件のもと逆解析を行ったところ，図-2に示す評価関数が得られた．評価関数は減少し，おおむね零へ収束していることがわかる．また，図-3および図-4に入力点における音圧波形の経時変化と，計測点における音圧波形の経時変化を示す．図-3により，逆解析の計算により得られた音圧波形（実線）は終端時刻付近において正解として設定した音圧波形（破線）と一致していない部分もあるが，それ以外の時刻においては，おおむね一致おり，正解値とほぼ等しい値が得られていることがわかる．
また，図-4は，計測点において得られた最終イタレーション時の計測点での音圧波形（実線）と正解として設定した音圧波形（破線）の比較を示している．

図-3

入力点における音圧波形の経時変化

図-4

計測点における音圧波形の経時変化

結果より，計測点においては，計算値と計測値はおおむね一致していることがわかる．音圧の計測点において，計算値（実線）と計測値（破線）が一致しているにも関わらず，入力点において計算値（実線）と正解値（破線）が終端時刻付近において一致しない要因としては，入力点における音圧が計測点に届くまでの時間差が影響しているものと考えられる．
このことより，逆解析により得られた入力点における音圧波形は，終端時刻付近においては，計測点における影響を含まずに求められていると考えられ，信頼性が薄い結果となっていると思われる．よって，本論文において紹介した逆解析により，入力点における音圧履歴を求める場合には，十分に計測値の計測時間を確保し，終端時刻付近の入力点における音圧変化は取り除いて，入力点における音圧履歴の同定値とすることが望ましいと考えられる．また，図-5～8に評価関数収束時における音圧の分布を示す．
100秒おきの結果を示しており，入力点において与えた音圧は適切に伝播し解析を行えていることを確認できる．

- 37 -

非破壊打音試験における打音信号同定のための逆解析

図-5

評価関数収束時における音圧の分布(100sec.)

図-8

評価関数収束時における音圧の分布(400sec.)

3.3 数値解析例２正弦波による逆解析ではおおむね一致した結果を得ることができた．しかし，実際の打音信号はパルス波に近く，その波形の大きさも打撃力や，対象とする材料により異なる．この様な場合，打音試験の際，予め打音信号を計測値から求めることが可能であれば，打音試験援用シミュレーションの際も利用できるものと考えられる．そこで本検討では，式(19)に示すパルス波による音圧波形を入力点における正解値として仮定し，正解となる音圧波形を逆解析により算定する検討を行う．式(19)において，Rはパルスの立ち上がる時間を示しており，R=50sec.
図-6

評価関数収束時におけると設定する．

音圧の分布(200sec.)p  exp  (t  R) 2



(19)

逆解析の結果，図-9に示す評価関数が得られた．前節と同様に，評価関数はおおむね零へ収束していることがわかる．また，図-10および図-11に入力点における音圧波形の経時変化と，計測点における音圧波形の経時変化を示す．図-10より，逆解析の計算により得られた音圧波形は，パルスの前後にリンギング（パルス近傍における振動）が見られることや，パルスの最大値が正解値では1.0であるのに対して0.42程度となっていることを確認できる．しかし，定性的には，パルスの立ち上がる時間はおおむね一致しており，正解値に近い値が得られているこ
図-7

評価関数収束時におけるとがわかる．また，図-11は，計測点において得ら音圧の分布(300sec.)

れた最終イタレーション時の計測点での音圧波形と- 38 -

倉橋 貴彦・松岡 瑛喜・村上 祐貴・外山 茂浩・池田 富士雄・井山 徹郎・田原 喜宏

正解として設定した音圧波形の比較を示している．
結果より，計測点においては，計算値と計測値はおおむね一致していることがわかる．また，図-12～15に評価関数収束時における音圧の分布を示す．
100秒おきの結果を示しており，入力点において与えた音圧は伝播し，適切に解析が行えていることを確認できる．

図-12 評価関数収束時における音圧の分布(100sec.)

図-9 評価関数の収束履歴

図-13 評価関数収束時における音圧の分布(200sec.)
図-10 入力点における音圧波形の経時変化

図-11 計測点における音圧波形の経時変化

図-14 評価関数収束時における音圧の分布(300sec.)- 39 -非破壊打音試験における打音信号同定のための逆解析

参考文献1)

浅野 雅則, 鎌田 敏郎, 国枝 稔, 六郷 恵哲, コンクリー
ト内部欠陥の寸法および深さと打音特性値の定量的関係 , コンク リート工学 年次論 文集 , Vol.23, No.1,pp.589-594, 2001.

2)

日本建築学会，音環境の数値シミュレーション－波
動音響解析の技法と応用－，丸善株式会社, 2011.

3)

S.Krenk, Dispersion-corrected explicit integration of thewave equation, Computer Methods in Applied Mechanics
and Engineering, Vol.191, pp.975-987, 2001.

4)

嘉納 秀明, システムの最適理論と最適化, コロナ社,
1987.

図-15 評価関数収束時における(2016. 8. 19 受付)

音圧の分布(400sec.)

４．結論
本論文では，非破壊打音試験における入力信号同定のための逆解析に関する検討を行った．音圧の伝播計算には波動方程式を適用し，空間方向の離散化には有限要素法，時間方向の離散化にはKrenkの提案する手法 3)を適用した．また，入力信号の同定問題に対しては随伴変数法を適用した．入力点において正解となる音圧波形を設定し，逆解析による数値実験を行った結果，同定された入力点における音圧波形は，終端時刻付近において，正解として仮定した音圧と異なる波形を示す結果となった．この結果は，入力点における音圧が計測点へ伝わる時間差が影響していると考えられる．そのため，本論文において紹介した入力点における音圧波形の同定を行う
場合は，十分な計測時間を確保し逆解析を行う必要があると考えられる．また，パルス波による逆解析においてもおおむね正しい結果を得ることができた．
将来的な検討としては，本解析プログラムを拡張し，配筋状態予測シミュレーションを行う予定である．

謝辞：本論文を執筆するにあたり「ダイハツ工業防振プロジェクト（長岡技術科学大学）」，「次世代ものづくり技術の基盤となる超高信頼性材料創成事業（長岡技術科学大学）」の援助を受けた．ここに謝意を表す．また，本論文の解析結果は，九州大学情報基盤研究開発センターの高性能アプリケーションサーバSR16000を使用し計算を行ったものである．
関係皆様へ謝意を表す．
- 40 -