JP3601907B2 - Vibration characteristic analyzer - Google Patents

Description

【００１１】
ここで、ｍは同時に参照するＦＲＦの数、Ｎは周波数点の数、そして、Ｗ_ｉは重み係数であり、Ｗ_ｉはｙ_ｉの母分散の逆数に比例する。
【００１２】
上記の残差二乗和Ｓが極小となる｛γ^＊｝＝｛γ｝＋｛Δγ｝をニュートン法で探索すれば次式で求めることができる。
【００１３】
【数２】

【００１４】
ここで［Ａ］は∂^２Ｓ／∂γ_ｉ∂γ_ｊを成分とする行列でヘシアン行列と呼ばれている。
【００１５】
従来の偏分反復法では、式（２）における行列［Ａ］と右辺を式（１）から求め、式（２）を解くことにより修正ベクトル｛Δγ｝を求める。そして｛γ｝の初期値から始めて、繰り返し修正ベクトル｛Δγ｝を求めることによりＳが極小値をとるときのモード特性を算出している。
【００１６】
この方法では、ＦＲＦの数ｍが増加すると｛γ｝の元数が増加するために、計算時間が急増する。また繰り返し計算において発散して解が得られない場合がある。したがって、従来の偏分反復法は主にＦＲＦを１個づつ個別に処理するのに使用されている。
【００１７】
そこで本発明者は、モード特性を求めるための上記のような周波数領域法の有する長所を保持したままで、次のような改良を提案した。
【００１８】
まず、特開平８−５４５０号公報において、非線形最小二乗法と最適化手法を合わせることにより初期値が真値と大きく異なるときの収束性と精度の大幅な改善を提案した。
【００１９】
すなわち、本発明者の手法では、求めるべきモード特性を固有モード形状｛α｝と固有振動数・モード減衰比｛β｝とに分ける。この場合、｛α｝はＦＲＦの理論式ｆの中に線形項として、また｛β｝は非線形項として含まれており、｛α｝を｛β｝の関数として表すことができる。これにより未知数は｛β｝のみになる。
【００２０】
このときＳが極小になる｛β^＊｝＝｛β｝＋｛Δβ｝は次式から求めら れる。
【００２１】
【数３】

【００２２】
ヘシアン行列［Ｂ］と右辺が算出できれば式（３）は計算できる。
【００２３】
また、Ｓの｛α｝に関する極小条件下で｛β｝を探索するために、収束性が改善される。｛β｝を探索するときに固有振動数とモード減衰比各々に適した側面制約を設定すれば、収束性はさらに改善される。
【００２４】
以上の方法により、実験で得られた多数のＦＲＦに共通で最も確からしい非線形項｛β｝を求めることができる。
【００２５】
非線形項｛β｝が求まれば、これを式（１）に代入し、残った線形項｛α｝を次式で求める。
【００２６】
【数４】

【００２７】
｛α｝は線形項なので式（４）は線形最小二乗法の式と一致し、繰り返し計算を必要としないで対象とする全周波数領域に対して｛α｝を求めることができる。
【００２８】
ＦＲＦの数をｍ、採用モード数をｎと表すと、式（４）における行列［Ｃ］の大きさは（２ｎ＋４）ｍである。
【００２９】
この場合、行列［Ｃ］の成分Ｃ_ｉｊ＝∂^２Ｓ／（∂α_ｉ∂α_ｊ）はｉとｊが異なるＦＲＦに属する場合には零となる。したがって、行列［Ｃ］は、ｍ個の（２ｎ＋４）×（２ｎ＋４）の行列に分解できるので計算時間と記憶容量が小さくて済む。
【００３０】
【発明が解決しようとする課題】
しかしながら、例えば、ディーゼルエンジンの低騒音化を目的にして、エンジンブロックの設計支援を実験モード解析で行うことを考える場合には、最低次の基本モードから騒音の原因となる高次モードまでの広い周波数領域に存在する多数のモード特性を実験で同定する必要がある。
【００３１】
近年のデータ処理装置の進歩により多数のスペクトル線を持つＦＲＦの収得が可能になり、現在６４０１本のスペクトル線を持つデータ処理装置が市販されている。このような装置でＦＲＦを騒音領域まで採取すると、これには多数の固有モードが含まれている。
【００３２】
この点に関しては、本発明者は特願平８−６８２５１号において、ＦＲＦの数を採用固有モード数にまで縮小することを提案した。このＦＲＦの数の縮小は採用固有モード数が２０個位までならば、計算の効率に大きく貢献する。
【００３３】
しかしながら、採用固有モード数がそれより増加すると縮小後のＦＲＦの数も増加し、式（３）の元数が増加するために、計算時間と記憶容量が急増するという問題点があった。
【００３４】
したがって、本発明は、被試験物を加振するためのインパルスハンマーに取り付けられた力検出用センサと、該被試験物の任意の場所に取り付けられて該被試験物の該インパルスハンマーの加振による応答を測定する測定センサと、該力検出用センサの出力信号及び該測定センサの出力信号とを受けて周波数応答関数データを求め、該周波数応答関数データと該周波数応答関数の理論式の残差二乗和を最小とする最小二乗法において線形項と非線形項に分けて該被試験物のモード特性を同定する演算装置と、を備えた振動特性解析装置において、対象とする全周波数領域に固有モードが多数存在する場合でも、計算時間と記憶容量を減少させることを目的とする。
【００３５】
【課題を解決するための手段】
〔１〕上記の目的を達成するために、本発明に係る振動特性解析装置は、演算装置が、周波数応答関数データにおける対象とする全周波数領域を、該周波数領域内の隣接する各固有モードのピーク周波数間の中間周波数を分割境界として固有モードが１個のみ存在する領域に分割し、さらに該隣接する領域間の分割境界に最も近い固有モード間の周波数と該隣接する領域内にそれぞれ存在する固有モード数との積に逆比例する値を引力と定義したときの隣接する領域間の全引力を計算し、該引力の最も大きな隣接した領域間を結合し、該結合された領域内の固有モード数が所定数を越えないように領域間結合を繰り返すことにより該分割された周波数領域の数を減少させ、該分割された周波数領域毎に該非線形項を演算し、該演算した非線形項により対象とする全周波数領域において線形項を演算することを特徴としている。
【００３６】
この同定法によれば、対象とする全周波数範囲を複数の対象領域に自動的に最適分割することにより、各分割領域内の固有モード数を同定に必要な計算時間と記憶容量が少なくて済むように規定した最大固有モード数以下にすることができる。これについて以下に詳しく説明する。
【００３７】
図１に示すように対象とする周波数領域を分割境界ａ及びｂで分割した場合は、領域Ｘ内に固有モードＡ，Ｂ，Ｃがあるとして計算した場合には、固有モードＣは領域Ｘ外でその近くに位置する固有モードＤの影響を受けて変形しているため計算値は不正確になる。
【００３８】
そこで、分割境界ｂの代わりに分割境界ｃで分割して領域Ｙ内に固有モードＡ〜Ｄがあるとして計算すれば領域Ｙ内の固有モードＢ，Ｃの相互の影響は考慮して計算されしかも領域Ｙから離れている固有モードＥによる影響は少なくて済むこととなる。
【００３９】
したがって、分割は隣接する固有モードのピーク周波数がなるべく離れている周波数点で行うことが同定誤差を少なくするために望ましい。また、分割後の各領域内の固有モード数を規定数以内でなるべく均等にして、分割領域内に多くの固有モードを入れて相互の影響を考慮した計算を実行して同定精度を上げることが望ましい。
【００４０】
そこで、本発明においては、まず分割された領域内に含まれる固有モードが１つだけになるように全周波数領域を細分割する。次に固有モードが最も近接した２個の隣接領域を結合することを繰り返す。結合後の領域内の固有モード数が規定数以上になる場合は領域間の結合が出来ないとする。
【００４１】
このような現象は自然界において水滴の凝集現象などと類似しており、結合後の形態も当初狙ったものに近いと考えられる。
【００４２】
以下に上記の原則を用いた最適分割の具体的な手順を図２のフローチャート及び図３のグラフを用いて説明する。
【００４３】
ステップＳ１：
対象とする周波数領域を各固有モードが１個になるように細分割する。具体的には、分割境界は、各固有モードのピーク周波数における隣接したもの同士間の中間周波数とする。
【００４４】
例えば、図３において、固有モードＡを有する周波数領域▲１▼と固有モードＢを有する周波数領域▲２▼とを分割するには、それらのピーク周波数ａ，ｂの中間周波数（ａ＋ｂ）／２を分割境界ｂ１とする。さらに、固有モードＢと固有モードＣとの間もそれらのピーク周波数の中間周波数（ｂ＋ｃ）／２を分割境界ｂ２とし、以下同様にして対象とする周波数領域は各々固有モードＡ〜Ｉを有する９個の周波数領域▲１▼〜▲９▼に分割される。
【００４５】
ステップＳ２：
次にステップＳ１で分割された領域間の引力Ｆを計算する。ここで、引力とは、隣接する固有モード間で相互に影響を与える要因の大きさを本発明において便宜的に定義したものである。引力Ｆは次の式（５）で与えられる。
【００４６】
【数５】
Ｆ＝１／（Ｌ×Ｎ１×Ｎ２） ・・・式（５）
Ｌ ：領域端の固有モード間の距離
Ｎ１：領域１内の固有モード数
Ｎ２：領域２内の固有モード数
【００４７】
図３を例にしてこの式（５）を説明すると、例えば固有モードＡを有する領域と固有モードＢを有する領域との間の引力Ｆ_ＡＢは、Ｆ_ＡＢ＝１／（Ｌ１×１×１）で与えられ、固有モードＢを有する領域と固有モードＣを有する領域との間の引力Ｆ_ＢＣは、Ｆ_ＢＣ＝１／（Ｌ２×１×１）で与えられる。
【００４８】
ステップＳ３：
次に、ステップＳ２で求めた各領域間の引力のうち最も大きい領域を結合する。
【００４９】
上記の例では、Ｆ_ＡＢが最も大きいものとすれば、固有モードＡを有する領域▲１▼と固有モードＢを有する領域▲２▼は２個の固有モードＡ，Ｂを有する領域Ｍとなる。
【００５０】
その結合される領域内の固有モード数が規定値（例えば８）を越える場合はその２つの領域は結合されない。
【００５１】
例えば、固有モードＡ，Ｂを有する領域Ｍと、固有モードＣ，Ｄを有する領域Ｎの間の引力Ｆ_ＭＮを求めると、領域端の固有モードＢ，Ｃ間の距離（周波数）はＬ_２で与えられ、領域Ｍには２個の固有モードＡ，Ｂと領域Ｎには２個の固有モードＣ，Ｄが含まれているのでＮ１，Ｎ２は共に２となり、Ｆ_ＭＮ＝１／（Ｌ_２×２×２）となる。
【００５２】
この引力Ｆ_ＭＮが他の領域間の引力よりも大きければ領域Ｍ，Ｎは結合されることになる。この場合には、この結合された領域の固有モード数は未だ４個であるから、規定数に達していないので結合することができる。
【００５３】
一方、領域Ｘと領域Ｙとが結合されないのは、領域Ｘの固有モード数が４個であり、領域Ｙの固有モード数が５個であるから、これらを結合すると、結合された領域に含まれる固有モード数「９」が上記の規定数「８」を越えるためである。
【００５４】
ステップＳ４：
ステップＳ２，Ｓ３の結合処理を繰り返して結合可能な領域が無くなった状態で終了する。
【００５５】
以上のようにして、周波数領域は自動分割される。
【００５６】
〔２〕また、上記の本発明〔１〕において、該演算装置が、該結合された領域内の固有モード数の規定数を該所定数より小さい数で該領域間結合を実行し、その後、該所定数で再度該領域間結合を実行することができる。
【００５７】
すなわち、固有モード数が大きい領域に挟まれた、固有モード数が１または２の領域が出現する可能性があり、この場合には対象とする周波数領域の分割数が減少せず該計算量の削減効果が小さくなってしまう。
【００５８】
これを避けるために、領域内の固有モードの規定数をｉと表すと、まずｉ／２を規定数として上記論理で結合し、その後にｉを規定数として結合してもよい。このようにすれば、固有モード数が１または２の領域が出現する可能性は大きく減少する。
【００５９】
〔３〕通常、実験で得られたＦＲＦには、対象領域外にも多数の固有モードが存在し、対象領域内のＦＲＦの振幅と位相に影響する。
【００６０】
この場合、本発明〔１〕の処理を実行した結果、図４（１）に示すような状態になったと仮定する。すなわち、分割境界ａ，ｂで分割された領域Ｘ内に３個の固有モードＡ，Ｂ，Ｃがあるとしてモード特性を同定すると、固有モードＡ，Ｂ，Ｃの山は、領域Ｘ外に固有モードが無いものとして同定されるので同図（２）において実線▲１▼で示す実際の形状より、破線▲２▼で示す如く固有モードＤの山の影響を受けて歪んでいる、この歪みは上記の本発明〔１〕での領域分割の仕方によっては顕著なものとなる。
【００６１】
このような誤差を補正するために、従来では剰余剛性Ｚと剰余質量Ｓという２定数を導入していた。これらの定数Ｚ，Ｓは、対象領域の境界から十分離れた固有モード群を表す係数として近似的に定義される。これらを用いた粘性減衰の場合のＦＲＦを次式に示す。
【００６２】
【数６】

【００６３】
ここで、Ｇ（ω）はコンプライアンスＦＲＦ、ωは角振動数、ｊは虚数単位、λｒ、Ａｒはｒ次固有モードの複素固有値と留数、ｎは採用固有モード数、＊は複素共役を表す。
【００６４】
試行計算の結果、上記の自動分割によるモード特性の同定の場合には、式（６）のように従来の剰余質量Ｓと剰余剛性Ｚを導入するだけではこの部分が誤差として影響して来るため精度的に不十分であることが判明した。
【００６５】
これを改良する手法としては、（ｉ）ローラン級数を利用する手法、及び（ｉｉ）領域外に固有モードを追加してこれで領域外の固有モード群を代表させる手法がある。
【００６６】
後者の手法（ｉｉ） は採用固有モード数＋２（領域の両外側）の自由度で系を表すことになり、前者の手法（ｉｉ）より意図が明快である。
【００６７】
上記の手法（ｉｉ）について図５を用いて説明する。
分割領域▲２▼においては、その領域両側において、実線で示すような振幅を本来有している筈であるが、上記の剰余質量と剰余剛性を用いただけでは、この部分が誤差として影響して来る。
【００６８】
そこで、この誤差を無くして精度のよい固有モード値を得るために、破線で示す擬似固有モードを上記の剰余質量と剰余剛性を含む上記のモード特性の同定式（６）に与えてやれば、これを補正することができる。
【００６９】
すなわち、手法（ｉｉ）と従来の剰余質量・剰余剛性を併用し、領域外の固有モードのうち、十分遠方の固有モード群を剰余質量と剰余剛性で表わし、また境界両側近傍の固有モード群を各々一つの擬似固有モードで表すものである。
【００７０】
この結果、剰余質量と剰余剛性を考慮した特性モードの同定式（６）に、さらに２つの擬似固有モードを加えると次式のようになる。
【００７１】
【数７】

【００７２】
ただし、領域外の固有モードの影響を補正するために、領域外に擬似固有モードを追加すると、繰り返し計算において収束性が悪化する。
【００７３】
この問題を回避するために、本発明ではモード特性の同定を次のように２段階に分ける。
【００７４】
すなわち、第１段階では、従来の剰余質量と剰余剛性だけを採用してモード特性を同定する。この同定方法は、収束性が良く、良好な精度の固有振動数が求まる。しかし、モード減衰は領域外の固有モードの影響を受けて誤差を含んでいることが多い。
【００７５】
第２段階では、領域外の固有モードの影響を補正するために領域外に２個の擬似固有モードを追加し、第１段階の結果を初期値にして固有振動数とモード減衰とを求める。この場合、固有振動数の初期値が第１段階で精度良く求まっているので、第２段階における繰り返し計算の収束は容易である。
【００７６】
そして、この補正演算を分割領域毎に行って目的とする固有モードが得られる。
【００７７】
【実施例】
図６は５１個の質点と５０個のバネ（粘性減衰素子を含む）の直列結合からなる力学モデルを示しており、３０００Ｈｚまでに１個の剛体モードと５０個の固有モードを有している。図７はこの力学モデルのＦＲＦ（アクセレランス）をグラフで示している。この力学モデルでは粘性減衰の仮定を採用し、そのモード減衰比は全固有モードで０．００２とする。
【００７８】
図７から分かるように周波数が高くなればなるほど、減衰力が大きくなると同時に固有振動数が近接するので、各固有モードが重なり合い、モード特性の同定が困難になることが予想される。
【００７９】
同定の対象周波数領域は５０〜２９００Ｈｚで、０．５Ｈｚ間隔の５７０１個の周波数点を採用している。また、本発明により自動的に最適分割された境界を破線で示す。分割領域内の固有モード数は一例として「８」以内と規定した。
【００８０】
その結果、対象領域は自動的に７個に分割されている。分割後の領域内の固有モードの数は低周波数側からそれぞれ６，４，８，８，８，８，８であり、規定数「８」以内に収まってる。分割後の各領域を低周波側から、領域▲１▼〜▲７▼と表す。
【００８１】
モード特性パラメータの一つである固有振動数の同定結果は、ここでは示さないが正解値に極めて近いものであった。
【００８２】
また、同じくモード特性パラメータの一つであるモード減衰比の同定結果を図８に示す。
【００８３】
×印は従来の剰余質量と剰余剛性のみで補正した結果（式（６））を表し、各分割領域において境界付近と中央付近で同定誤差が大きく同定精度が低下していることを示している。領域▲１▼の低周波側でやや精度が良いのは、このモデルが剛体モードを有するために、剰余質量の仮定が成立していることによる。
【００８４】
＊印は剰余質量と剰余剛性による補正を行わずに、領域外の低周波側と高周波側に各１個の擬似固有モードを挿入して補正した結果を示す。この場合には、１５００Ｈｚ以下の低周波数領域において同定誤差が小さく良好な結果を得ているが、１５００Ｈｚ以上の高周波数領域では同定誤差が増加していることが分かる。
【００８５】
これは高周波数領域では図７に示すように固有モード間の間隔が狭くなり、分割境界の外のすぐ近くに隣接領域の固有モードが存在し、この近接した固有モードの影響と遠方の固有モードの影響を合わせて一つの固有モードで補正することが困難になるためと考えられる。
【００８６】
○印は本発明による手法で、剰余質量と剰余剛性による補正に、領域の両側にそれぞれ１個の擬似固有モードを追加した例を表す。この場合には、領域▲６▼まで良好な同定結果である。領域▲７▼の高周波側で精度が低下しているのは、図７から分かるように固有モードが極めて近接していて分離が困難なためである。
【００８７】
図９には固有モード（図６の力学モデルの左端の質点の振幅）の同定誤差を示す。この例では図７と同様の傾向を示す。ここでも○印で示す本発明の手法は良好な同定結果を示している。
【００８８】
図１０と図１１に計算時間と記憶容量の比較グラフが示されている。×印は剰余質量と剰余剛性のみによる補正で、本発明による自動最適分割を行わない場合を示し、○印は本発明による手法で補正し対象領域を自動最適分割した場合を示す。
【００８９】
これより、固有モード数が１５以上になると、対象領域を自動最適分割しない手法（×印）では計算時間と記憶容量が共に急増し、採用固有モード数４５以上では計算出来なかったことが分かる。
【００９０】
また、固有モード数４０で比較すると、対象領域を自動最適分割することにより、計算時間は１／５（図１０）に、記憶容量は１／９（図１１）に減少していることが分かる。
【００９１】
採用固有モードが５０以上に増加しても、自動最適分割を行えば、計算時間の増加は緩やかであり、記憶容量はほどんど増加していない。
【００９２】
このように、本発明手法によれば、多数のモード特性を、実用的な計算時間と記憶容量で精度良く求めることができる。
【００９３】
なお、上記の実施例において、対象周波数領域内における固有モードを同定する場合、式（１）で示す残差二乗和Ｓは固有振動数がわずかに変化しても大きく変化する。そこで、特開平８−５４５０号公報においては、繰り返し計算において、固有振動数にはその変化率が０．０５以下、モード減衰比には０．８以下という側面制約を設定し、発散減少を防止している。
【００９４】
従って本発明において、領域外の２個の擬似固有モードに対しては、固有振動数とモード減衰比の変化率が共に０．１以下という側面制約を設定して収束の遅れを防いでもよい。
【００９５】
さらに、領域外に擬似固有モードを追加した場合、式（３）の行列［Ｂ］を算出するために式（４）の行列［Ｃ］の逆行列を使用すると、［Ｃ］の特異解分解の結果、この逆行列の計算がランク落ちにより不可能になる可能性がある。
【００９６】
従って、逆行列を用いないガウスの消去法による直接解法に改良することも可能である。
【００９７】
【発明の効果】
以上説明したように、本発明に係る振動特性解析装置によれば、演算装置が、周波数応答関数データにおける対象とする全周波数領域を、固有モード数が所定数を越えない周波数領域に自動最適分割し、該分割された周波数領域毎に該非線形項を演算し、該演算した非線形項により対象とする全周波数領域において線形項を演算することにより、多数の固有モードのモード特性を同定する場合に、計算時間と記憶容量が大きく低減することが可能となった。
【００９８】
また、該演算装置が、分割された各周波数領域外の遠方の固有モードを剰余質量及び剰余剛性を用いて表すことによりモード特性を演算し、その後に該演算結果を初期値として、該分割された周波数領域の両側の周波数領域に各々補正用の擬似固有モードを与えて該非線形項及び該線形項を演算することにより、対象領域外の固有モードの影響を、精度良く補正することが可能となった
【図面の簡単な説明】
【図１】本発明に係る振動特性解析装置の演算装置による周波数領域分割例を示したグラフ図である。
【図２】本発明に係る振動特性解析装置の演算装置による自動最適周波数領域分割を説明するためのフローチャート図である。
【図３】本発明に係る振動特性解析装置の演算装置による自動最適周波数領域分割例を示したグラフ図である。
【図４】対象周波数領域外の固有モードが対象周波数領域内の固有モードに与える影響を示したグラフ図である。
【図５】本発明に係る振動特性解析装置の演算装置による対象周波数領域外の固有モードによる補正手法を説明するためのグラフ図である。
【図６】本発明に係る振動特性解析装置の演算装置における有効性の検証で用いる力学モデルを示したブロック図である。
【図７】図６の力学モデルにおけるＦＲＦ（アクセレランス）の振幅を示すグラフ図である。
【図８】図６の力学モデルにおけるモード減衰比の同定誤差を示すグラフ図である。
【図９】図６の力学モデルにおける固有モード振幅の同定誤差を示すグラフ図である。
【図１０】図６の力学モデルにおけるモード特性同定に要する計算時間を示すグラフ図である。
【図１１】図６の力学モデルにおけるモード特性同定に必要とする記憶容量を示すグラフ図である。
【図１２】振動特性解析装置の一般的な構成例を示したブロック図である。
【符号の説明】
１ 被試験物
２ インパルスハンマー
３ 力検出用センサ
４ 測定センサ
５ 演算装置
図中、同一符号は同一又は相当部分を示す。TECHNICAL FIELD OF THE INVENTION
The present invention relates to a vibration characteristic analyzing apparatus, and more particularly to a vibration characteristic analyzing apparatus that obtains a response by exciting a test object such as an engine block and identifies a mode characteristic of the test object from a response function.
[0001]
In order to analyze the dynamic characteristics of an actual mechanical structure, it is necessary to identify the dynamic characteristics of the structure from the vibration response characteristics obtained by performing a vibration test on the structure. At present, the mainstream method is to identify modal parameters (modal parameters such as natural frequency, mode damping ratio, and eigenmode shape), which are parameters that directly express the phenomenon, are easy to understand the phenomenon directly, and have good convergence. ing.
[0002]
[Prior art]
FIG. 12 schematically shows a configuration of a conventionally known vibration characteristic analyzing apparatus. First, a DUT 1 and an impulse hammer 2 for exciting the DUT 1 are prepared.
[0003]
Then, a force detecting sensor 3 is attached to the hammer 2, and a three-dimensional (x, y, z) direction acceleration of a vibration response due to the vibration of the hammer 2 is detected at an arbitrary position of the test object 1. Is attached.
[0004]
Then, an output signal of the acceleration sensor 4 and an output signal of the force detection sensor 3 are provided to the arithmetic unit 5 to perform FFT processing, obtain frequency response function data (FRF: compliance), calculate mode characteristics, and output the result. are doing.
[0005]
In this case, the vibration experiment is usually performed while the vibration place of the DUT 1 by the hammer 2 is fixed, and the acceleration sensor 4 is sequentially moved to give a plurality of response functions to the arithmetic unit 5. The DUT 1 is fixed by a soft spring (not shown) or the like.
[0006]
Many algorithms have been proposed as algorithms for the arithmetic unit (experimental mode analyzer) 5 of such a vibration characteristic analyzer. Among them, the "partial repetition method" in the frequency domain has high accuracy. Currently widely used.
[0007]
This method identifies modal characteristics by performing iterative calculation while testing an error, and can obtain accurate modal characteristics.
[0008]
However, the calculation efficiency is low, and a large amount of calculation time is required particularly when many frequency response functions are referred to at the same time. If the initial value is not properly given, the calculation may diverge and a solution may not be obtained.
[0009]
That is, when the mode characteristics {gamma}, in order to obtain the {gamma}, the experimental data of the FRF at the frequency point (ω _i) y _i, the value of the corresponding theoretical formula When expressed by f _i, mixed If the resulting noise is due to a number of independent random factors, the distribution of y _i will be a normal distribution. In this case, the most likely mode characteristic is {γ} that minimizes the weighted residual sum of squares S of the following equation.
[0010]
(Equation 1)

[0011]
Here, m is the number of FRF referring simultaneously, N is the number of frequency points and,, W _i is a weighting factor, W _i is proportional to the inverse of the population variance of y _i.
[0012]
If {γ ^* } = {γ} + {Δγ} at which the above-mentioned residual sum of squares S is minimized is found by the Newton method, it can be obtained by the following equation.
[0013]
(Equation 2)

[0014]
Here, [A] is a matrix having ∂ ² S / ∂γ _i ∂γ _j as components and is called a Hessian matrix.
[0015]
In the conventional partial repetition method, the matrix [A] and the right side in Expression (2) are obtained from Expression (1), and the correction vector {Δγ} is obtained by solving Expression (2). Then, starting from the initial value of {γ}, the mode characteristic when S takes the minimum value is calculated by repeatedly calculating the correction vector {Δγ}.
[0016]
In this method, when the number m of FRFs increases, the element number of {γ} increases, so that the calculation time rapidly increases. In some cases, a solution cannot be obtained due to divergence in the iterative calculation. Therefore, the conventional iterative iteration method is mainly used to process FRFs individually one by one.
[0017]
Therefore, the present inventor has proposed the following improvements while retaining the advantages of the above-described frequency domain method for obtaining the mode characteristics.
[0018]
First, Japanese Patent Laid-Open Publication No. Hei 8-5450 proposes a significant improvement in convergence and accuracy when the initial value is significantly different from the true value by combining the nonlinear least squares method and the optimization method.
[0019]
That is, in the method of the inventor, the mode characteristic to be obtained is divided into the eigenmode shape {α} and the eigenfrequency / mode damping ratio {β}. In this case, {α} is included as a linear term in the theoretical formula f of the FRF, and {β} is included as a nonlinear term, so that {α} can be expressed as a function of {β}. Thus, the only unknown is {β}.
[0020]
At this time, {β ^* } = {β} + {Δβ} at which S is minimized is obtained from the following equation.
[0021]
(Equation 3)

[0022]
If the Hessian matrix [B] and the right side can be calculated, equation (3) can be calculated.
[0023]
Further, since {β} is searched under the minimum condition regarding {α} of S, the convergence is improved. Convergence is further improved by setting side constraints suitable for the natural frequency and the mode damping ratio when searching for {β}.
[0024]
By the above method, the most probable nonlinear term {β} common to many FRFs obtained in the experiment can be obtained.
[0025]
Once the nonlinear term {β} is determined, it is substituted into equation (1), and the remaining linear term {α} is determined by the following equation.
[0026]
(Equation 4)

[0027]
Since {α} is a linear term, equation (4) matches the equation of the linear least squares method, and {α} can be obtained for the entire frequency range of interest without the need for repeated calculations.
[0028]
Assuming that the number of FRFs is m and the number of adopted modes is n, the size of the matrix [C] in equation (4) is (2n + 4) m.
[0029]
In this case, the component C _ij = ∂ ² S / (∂α _i ∂α _j ) of the matrix [C] becomes zero when i and j belong to different FRFs. Therefore, since the matrix [C] can be decomposed into m (2n + 4) × (2n + 4) matrices, the calculation time and storage capacity can be reduced.
[0030]
[Problems to be solved by the invention]
However, for example, in order to reduce the noise of a diesel engine and to consider design support for an engine block by an experimental mode analysis, a wide range from the lowest-order basic mode to a higher-order mode that causes noise is considered. It is necessary to identify a large number of mode characteristics existing in the frequency domain through experiments.
[0031]
Recent advances in data processing devices have made it possible to obtain FRFs having a large number of spectral lines, and data processing devices having 6401 spectral lines are currently on the market. When the FRF is collected to the noise region by such a device, it includes a large number of eigenmodes.
[0032]
In this regard, the present inventor has proposed in Japanese Patent Application No. Hei 8-68251 to reduce the number of FRFs to the number of adopted eigenmodes. This reduction in the number of FRFs greatly contributes to the efficiency of calculation if the number of adopted eigenmodes is up to about 20.
[0033]
However, when the number of adopted eigenmodes increases, the number of FRFs after reduction also increases, and the number of elements in equation (3) increases, so that there is a problem that the calculation time and the storage capacity rapidly increase.
[0034]
Therefore, the present invention provides a force detection sensor attached to an impulse hammer for exciting a test object, and an excitation of the impulse hammer for the test object attached to an arbitrary position of the test object. A frequency response function data by receiving the output signal of the force detection sensor and the output signal of the measurement sensor, and obtains the frequency response function data and the remainder of the theoretical equation of the frequency response function. An arithmetic unit that identifies a modal characteristic of the device under test by dividing into a linear term and a non-linear term in a least square method that minimizes the sum of squared differences. It is an object to reduce the calculation time and the storage capacity even when there are many modes.
[0035]
[Means for Solving the Problems]
[1] In order to achieve the above object, in a vibration characteristic analyzing apparatus according to the present invention, an arithmetic unit calculates all target frequency regions in frequency response function data for each adjacent eigenmode in the frequency region. The intermediate frequency between the peak frequencies is used as a division boundary to divide the region into which only one eigenmode exists, and further, the frequency between the eigenmodes closest to the division boundary between the adjacent regions and the frequency existing in the adjacent region When the value inversely proportional to the product of the number of eigenmodes is defined as the attractive force, the total attractive force between the adjacent regions is calculated, the adjacent regions having the largest attractive force are connected, and the eigenvalue in the connected region is calculated. By repeating the inter-region coupling so that the number of modes does not exceed a predetermined number, the number of the divided frequency regions is reduced, and the nonlinear term is calculated for each of the divided frequency regions. It is characterized by computing the linear term in the entire frequency range of interest by.
[0036]
According to this identification method, the calculation time and storage capacity required for identifying the number of eigenmodes in each divided region can be reduced by automatically optimally dividing the entire target frequency range into a plurality of target regions. Thus, the number of eigenmodes can be made equal to or less than the maximum number of eigenmodes defined as above. This will be described in detail below.
[0037]
As shown in FIG. 1, when the target frequency domain is divided by the division boundaries a and b, when the calculation is performed assuming that the eigenmodes A, B, and C are present in the area X, the eigenmode C is outside the area X. , Is deformed under the influence of the eigenmode D located in the vicinity, and the calculated value becomes inaccurate.
[0038]
Therefore, if the calculation is performed on the assumption that there are eigenmodes A to D in the area Y by dividing the area by the division boundary c instead of the division boundary b, the calculation is performed in consideration of the mutual influence of the eigenmodes B and C in the area Y. The influence of the eigenmode E which is far from the area Y is small.
[0039]
Therefore, it is desirable to perform the division at frequency points where the peak frequencies of adjacent eigenmodes are as far apart as possible to reduce identification errors. It is also possible to increase the identification accuracy by making the number of eigenmodes within each of the divided regions equal to or less than a specified number as much as possible, and performing calculations in consideration of mutual influences by putting many eigenmodes in the divided regions. desirable.
[0040]
Therefore, in the present invention, first, the entire frequency region is subdivided so that only one eigenmode is included in the divided region. Next, the coupling of the two adjacent regions having the closest eigenmode is repeated. If the number of eigenmodes in the combined area is equal to or greater than a prescribed number, it is assumed that the connection between the areas cannot be made.
[0041]
Such a phenomenon is similar to the aggregation phenomenon of water droplets in the natural world, and the morphology after bonding is considered to be close to the target initially.
[0042]
Hereinafter, a specific procedure of the optimal division using the above principle will be described with reference to the flowchart of FIG. 2 and the graph of FIG.
[0043]
Step S1:
The target frequency domain is subdivided so that each eigenmode becomes one. Specifically, the division boundary is an intermediate frequency between adjacent ones at the peak frequency of each eigenmode.
[0044]
For example, in FIG. 3, in order to divide the frequency region (1) having the eigenmode A and the frequency region (2) having the eigenmode B, the intermediate frequency (a + b) / 2 of the peak frequencies a and b is calculated. It is assumed that the division boundary is b1. Further, also between the eigenmode B and the eigenmode C, the intermediate frequency (b + c) / 2 of their peak frequencies is set as a division boundary b2, and similarly, the target frequency regions have the eigenmodes A to I, respectively. The frequency domain is divided into (1) to (9).
[0045]
Step S2:
Next, an attractive force F between the divided regions is calculated in step S1. Here, the gravitational force is a value defined for the sake of convenience in the present invention for the size of a factor that mutually affects adjacent eigenmodes. The attractive force F is given by the following equation (5).
[0046]
(Equation 5)
F = 1 / (L × N1 × N2) Equation (5)
L: distance between eigenmodes at the end of the region N1: number of eigenmodes in region 1 N2: number of eigenmodes in region 2
This equation (5) will be described with reference to FIG. 3 as an example. For example, the attractive force F _AB between the region having the eigenmode A and the region having the eigenmode B is F _AB = 1 / (L1 × 1 × 1) given attractive force _{F BC} between the region having an area and eigenmodes C with eigenmode B _is given by F BC = 1 / (L2 × 1 × 1).
[0048]
Step S3:
Next, the largest area of the attraction between the areas obtained in step S2 is combined.
[0049]
In the above example, _assuming that F _AB is the largest, the area (1) having the eigenmode A and the area (2) having the eigenmode B become an area M having two eigenmodes A and B.
[0050]
If the number of eigenmodes in the combined area exceeds a specified value (for example, 8), the two areas are not combined.
[0051]
For example, a region M having a eigenmode A, B, eigenmodes C, and obtaining the attraction F _MN between regions N with D, eigenmodes B region end, the distance between C (frequency) by L ₂ given area M 2 pieces of eigenmodes a in, B and area N 2 pieces of eigenmodes C in, because it contains D N1, N2 are both 2 _{becomes, F MN = 1 / (L} 2 × 2 × 2).
[0052]
If the attractive force _FMN is larger than the attractive force between the other regions, the regions M and N are combined. In this case, since the number of eigenmodes in the combined area is still four, it does not reach the specified number, so that the area can be combined.
[0053]
On the other hand, the reason that the region X and the region Y are not combined is that the region X has four eigenmodes and the region Y has five eigenmodes. This is because the number of eigenmodes “9” exceeds the specified number “8”.
[0054]
Step S4:
The combining process of steps S2 and S3 is repeated, and the process ends when there is no more combineable area.
[0055]
As described above, the frequency domain is automatically divided.
[0056]
[2] Also, in the above present invention [1], the arithmetic unit executes the inter-region coupling with a prescribed number of eigenmodes in the coupled region smaller than the predetermined number, and thereafter, The inter-region coupling can be executed again with the predetermined number.
[0057]
In other words, there is a possibility that a region having the number of eigenmodes of 1 or 2 sandwiched between regions having a large number of eigenmodes may occur. In this case, the number of divisions of the target frequency region does not decrease and the amount of calculation does not increase. The reduction effect is small.
[0058]
In order to avoid this, if the prescribed number of the eigenmodes in the area is represented by i, first, i / 2 may be used as the prescribed number and the above-described logic may be used for coupling, and then i may be used as the prescribed number. In this way, the possibility that an area having the number of eigenmodes of 1 or 2 appears is greatly reduced.
[0059]
[3] Normally, in an FRF obtained by an experiment, a large number of eigenmodes exist outside the target region, which affects the amplitude and phase of the FRF in the target region.
[0060]
In this case, it is assumed that the state shown in FIG. 4A is obtained as a result of executing the processing of the present invention [1]. That is, when the mode characteristics are identified assuming that there are three eigenmodes A, B, and C in the area X divided by the division boundaries a and b, the peaks of the eigenmodes A, B, and C are outside the area X. Since it is identified as having no mode, it is distorted by the influence of the peak of the eigenmode D as shown by the broken line (2) from the actual shape shown by the solid line (1) in FIG. This becomes remarkable depending on the manner of area division in the present invention [1].
[0061]
In order to correct such an error, two constants, that is, a residual stiffness Z and a residual mass S, have been conventionally introduced. These constants Z and S are approximately defined as coefficients representing a group of eigenmodes sufficiently distant from the boundary of the target area. The FRF in the case of viscous damping using these is shown in the following equation.
[0062]
(Equation 6)

[0063]
Here, G (ω) is the compliance FRF, ω is the angular frequency, j is the imaginary unit, λr, Ar is the complex eigenvalue and residue of the rth-order eigenmode, n is the number of adopted eigenmodes, and * is the complex conjugate. .
[0064]
As a result of the trial calculation, in the case of the identification of the mode characteristics by the above-mentioned automatic division, simply introducing the conventional residual mass S and the residual stiffness Z as shown in the equation (6) affects this part as an error. It turned out to be insufficient in accuracy.
[0065]
As a method for improving this, there are (i) a method using a Laurent series, and (ii) a method of adding an eigenmode outside the region and representing the eigenmode group outside the region.
[0066]
The latter method (ii) expresses the system with the degree of freedom of the number of adopted eigenmodes + 2 (both outside of the region), and the intention is clearer than the former method (ii).
[0067]
The above method (ii) will be described with reference to FIG.
In the divided region (2), both sides of the region should originally have an amplitude as shown by a solid line. However, if only the above-mentioned surplus mass and surplus stiffness are used, this portion will affect as an error. come.
[0068]
Therefore, in order to eliminate the error and obtain an accurate eigenmode value, a pseudo eigenmode indicated by a broken line is given to the above-described mode characteristic identification equation (6) including the residual mass and the residual stiffness. This can be corrected.
[0069]
That is, by using the method (ii) and the conventional residual mass / remainder stiffness together, the eigenmodes that are sufficiently far out of the eigenmodes outside the region are represented by the residual mass and the residual stiffness, and the eigenmodes near both sides of the boundary are represented. Each is represented by one pseudo eigenmode.
[0070]
As a result, the following expression is obtained by adding two pseudo eigenmodes to the identification expression (6) of the characteristic mode in consideration of the residual mass and the residual rigidity.
[0071]
(Equation 7)

[0072]
However, if a pseudo eigenmode is added outside the region in order to correct the influence of the eigenmode outside the region, the convergence deteriorates in the repeated calculation.
[0073]
In order to avoid this problem, in the present invention, the identification of the mode characteristic is divided into two stages as follows.
[0074]
That is, in the first stage, the mode characteristics are identified using only the conventional surplus mass and surplus rigidity. This identification method has good convergence and obtains a natural frequency with good accuracy. However, the mode attenuation often includes an error due to the influence of the eigenmode outside the region.
[0075]
In the second stage, two pseudo eigenmodes are added outside the region to correct the effects of the eigenmodes outside the region, and the natural frequency and mode damping are obtained using the results of the first stage as initial values. In this case, since the initial value of the natural frequency is accurately determined in the first stage, the convergence of the repeated calculation in the second stage is easy.
[0076]
Then, this correction operation is performed for each of the divided areas to obtain a desired eigenmode.
[0077]
【Example】
FIG. 6 shows a dynamic model composed of a series connection of 51 mass points and 50 springs (including viscous damping elements), and has one rigid mode and 50 eigenmodes up to 3000 Hz. . FIG. 7 is a graph showing the FRF (acceleration) of this dynamic model. This dynamic model adopts the assumption of viscous damping, and its mode damping ratio is set to 0.002 for all eigenmodes.
[0078]
As can be seen from FIG. 7, the higher the frequency, the greater the damping force and the closer the natural frequency, so that each natural mode is expected to overlap, making it difficult to identify the mode characteristics.
[0079]
The frequency range to be identified is 50 to 2900 Hz, and 5701 frequency points at 0.5 Hz intervals are employed. Also, a boundary automatically optimally divided by the present invention is indicated by a broken line. As an example, the number of eigenmodes in the divided area is defined to be within “8”.
[0080]
As a result, the target area is automatically divided into seven. The number of eigenmodes in the area after division is 6, 4, 8, 8, 8, 8, 8 from the low frequency side, respectively, and falls within the specified number “8”. Each area after division is represented as areas {circle around (1)} to {circle around (7)} from the low frequency side.
[0081]
Although the identification result of the natural frequency, which is one of the mode characteristic parameters, is not shown here, it was extremely close to the correct answer value.
[0082]
FIG. 8 shows the identification result of the mode attenuation ratio, which is also one of the mode characteristic parameters.
[0083]
The crosses indicate the result (Equation (6)) corrected only by the conventional residual mass and residual rigidity, and indicate that the identification error is large near the boundary and near the center in each divided region, and the identification accuracy is reduced. . The reason why the accuracy is slightly higher on the low frequency side of the region (1) is that the assumption of the residual mass is established because the model has a rigid body mode.
[0084]
The * mark shows the result of correction by inserting one pseudo eigenmode on each of the low frequency side and the high frequency side outside the region without performing the correction based on the residual mass and the residual rigidity. In this case, the identification error is small and a good result is obtained in a low frequency region of 1500 Hz or less, but the identification error increases in a high frequency region of 1500 Hz or more.
[0085]
This is because, as shown in FIG. 7, in the high frequency region, the interval between the eigenmodes becomes narrow, and the eigenmode of the adjacent region exists immediately outside the division boundary, and the influence of this eigenmode and the distant eigenmode It is considered that it is difficult to perform correction in one eigenmode in consideration of the influence of the above.
[0086]
The circles represent examples in which one pseudo eigenmode is added to each side of the region in addition to the correction based on the residual mass and the residual rigidity in the method according to the present invention. In this case, a good identification result is obtained up to the area (6). The accuracy is reduced on the high frequency side of the area (7) because the eigenmodes are extremely close to each other and difficult to separate as can be seen from FIG.
[0087]
FIG. 9 shows the identification error of the eigenmode (the amplitude of the leftmost mass point of the dynamic model of FIG. 6). This example shows the same tendency as FIG. Again, the method of the present invention indicated by a circle shows a good identification result.
[0088]
FIGS. 10 and 11 show comparison graphs of calculation time and storage capacity. A mark “X” indicates a case where the correction is performed only by the residual mass and the residual rigidity, and the automatic optimal division according to the present invention is not performed, and a mark “○” indicates a case where the correction is performed by the method according to the present invention and the target area is automatically optimally divided.
[0089]
From this, it can be seen that when the number of eigenmodes becomes 15 or more, both the calculation time and the storage capacity increase rapidly with the method of not automatically dividing the target region (marked by x), and the calculation cannot be performed with the number of adopted eigenmodes of 45 or more.
[0090]
Comparing with the number of eigenmodes of 40, it can be seen that the calculation time has been reduced to 1/5 (FIG. 10) and the storage capacity has been reduced to 1/9 (FIG. 11) by automatically and optimally dividing the target area. .
[0091]
Even if the adoption eigenmode increases to 50 or more, if the automatic optimal division is performed, the increase in the calculation time is gradual, and the storage capacity is hardly increased.
[0092]
As described above, according to the method of the present invention, a large number of mode characteristics can be accurately obtained with a practical calculation time and storage capacity.
[0093]
In the above embodiment, when the eigenmode in the target frequency region is identified, the residual sum of squares S represented by the equation (1) greatly changes even if the natural frequency slightly changes. Therefore, in Japanese Unexamined Patent Publication No. Hei 8-5450, in the repetitive calculation, the natural frequency is set to have a change rate of 0.05 or less and the mode damping ratio is set to 0.8 or less to prevent the divergence from decreasing. are doing.
[0094]
Therefore, in the present invention, for two pseudo eigenmodes outside the region, a side constraint that both the natural frequency and the rate of change of the mode damping ratio are 0.1 or less may be set to prevent delay of convergence.
[0095]
Further, when a pseudo eigenmode is added outside the region, when the inverse matrix of the matrix [C] of the equation (4) is used to calculate the matrix [B] of the equation (3), the singular solution decomposition of [C] is performed. As a result, the calculation of the inverse matrix may not be possible due to rank drop.
[0096]
Therefore, it is possible to improve the direct solution by Gaussian elimination without using the inverse matrix.
[0097]
【The invention's effect】
As described above, according to the vibration characteristic analysis apparatus of the present invention, the arithmetic unit automatically divides the entire target frequency region in the frequency response function data into frequency regions in which the number of eigenmodes does not exceed a predetermined number. When calculating the non-linear terms for each of the divided frequency domains, and calculating the linear terms in the entire frequency domain of interest by the calculated non-linear terms, when identifying the mode characteristics of a large number of eigenmodes, Thus, the calculation time and the storage capacity can be greatly reduced.
[0098]
Further, the arithmetic unit calculates the mode characteristic by expressing the distant eigenmodes outside the divided frequency regions using the residual mass and the residual stiffness, and thereafter, using the calculation result as an initial value, By applying pseudo-eigenmodes for correction to the frequency domains on both sides of the frequency domain and calculating the nonlinear term and the linear term, it is possible to accurately correct the influence of the eigenmodes outside the target area. No longer [Brief description of drawings]
FIG. 1 is a graph showing an example of frequency domain division by an arithmetic unit of a vibration characteristic analyzing apparatus according to the present invention.
FIG. 2 is a flowchart for explaining automatic optimal frequency domain division by a calculation device of the vibration characteristic analysis device according to the present invention.
FIG. 3 is a graph showing an example of automatic optimal frequency domain division by an arithmetic unit of the vibration characteristic analyzing apparatus according to the present invention.
FIG. 4 is a graph showing the effect of an eigenmode outside the target frequency domain on an eigenmode within the target frequency domain.
FIG. 5 is a graph for explaining a correction method based on an eigenmode outside a target frequency region by an arithmetic unit of the vibration characteristic analysis apparatus according to the present invention.
FIG. 6 is a block diagram showing a dynamic model used for verifying the validity of the arithmetic device of the vibration characteristic analysis device according to the present invention.
FIG. 7 is a graph showing the amplitude of FRF (acceleration) in the dynamic model of FIG. 6;
FIG. 8 is a graph showing an identification error of a mode damping ratio in the dynamic model of FIG. 6;
9 is a graph illustrating an identification error of an eigenmode amplitude in the dynamic model of FIG. 6;
FIG. 10 is a graph showing a calculation time required for identifying a mode characteristic in the dynamic model of FIG. 6;
11 is a graph showing a storage capacity required for identifying a mode characteristic in the dynamic model of FIG. 6;
FIG. 12 is a block diagram illustrating a general configuration example of a vibration characteristic analysis device.
[Explanation of symbols]
REFERENCE SIGNS LIST 1 test object 2 impulse hammer 3 force detection sensor 4 measurement sensor 5 arithmetic unit

Claims (3)

A force detection sensor attached to an impulse hammer for vibrating the DUT, and a measurement attached to an arbitrary position of the DUT to measure a response of the DUT by the excitation of the impulse hammer A sensor, receiving the output signal of the force detection sensor and the output signal of the measurement sensor, obtains frequency response function data, and minimizes the sum of squared residuals of the frequency response function data and the theoretical formula of the frequency response function. A computing device that identifies a mode characteristic of the DUT by dividing into a linear term and a non-linear term in a least-squares method,
The arithmetic unit sets the entire frequency region of interest in the frequency response function data to a region where only one eigenmode exists with an intermediate frequency between peak frequencies of adjacent eigenmodes in the frequency region as a division boundary. Divide, and further define the value which is inversely proportional to the product of the frequency between the eigenmodes closest to the division boundary between the adjacent regions and the number of eigenmodes present in the adjacent regions as the attractive force. Calculating the total attractive force between the regions, connecting the adjacent regions having the largest attractive force, and repeating the interregion coupling so that the number of eigenmodes in the coupled region does not exceed a predetermined number. Reducing the number of frequency domains obtained, calculating the non-linear terms for each of the divided frequency domains, and calculating the linear terms in all target frequency domains by the calculated non-linear terms. And the vibration characteristics analysis device. In claim 1,
The arithmetic unit executes the inter-region coupling with a prescribed number of the eigenmode numbers in the coupled region smaller than the predetermined number, and thereafter executes the inter-region coupling again with the predetermined number. Vibration characteristic analyzer that is a feature. In claim 1 or 2,
The arithmetic device calculates the mode characteristic by expressing the distant eigenmodes outside each of the divided frequency regions using the residual mass and the residual stiffness, and thereafter, using the calculation result as an initial value, A vibration characteristic analyzing apparatus, wherein a pseudo eigenmode for correction is given to each of frequency regions on both sides of the region to calculate the nonlinear term and the linear term.

Images