JP6178251B2

JP6178251B2 - ターボ機械の動特性演算方法およびターボ機械の動特性演算装置

Info

Publication number: JP6178251B2
Application number: JP2014013550A
Authority: JP
Inventors: 広敏在原; 弘毅佐成; 祥孝馬場; 裕樹亀山; 岡田　徹; 徹岡田
Original assignee: Kobe Steel Ltd
Current assignee: Kobe Steel Ltd
Priority date: 2013-10-25
Filing date: 2014-01-28
Publication date: 2017-08-09
Anticipated expiration: 2034-01-28
Also published as: JP2015108607A

Description

本発明は、ターボ機械のロータ系を構成する各機械要素の動特性を演算する技術に関する。

回転軸に連結された羽根車を回転させることで流体を圧縮または圧送するターボ機械において、高圧流体の流体不安定化力によりロータに自励振動が生じ、ターボ機械を停止させなければならない状況となることがある。ロータに自励振動が生じるか否かは、回転軸の剛性や質量、軸受の剛性や減衰係数、羽根車やシールで発生する不安定化力など、それぞれの機械要素の動特性を把握することによって適切に判断することができる。

ここで、特許文献１では、回転軸を加振することができるように構成された回転機械において、回転数を一定としつつ負荷条件を変化させ、各負荷条件ごとの振動応答データから回転機械の減衰比または許容励振係数を求めることで、回転機械の動特性を演算している。また、特許文献２にも、特許文献１と同様に、回転機械の減衰比または許容励振係数を求めることで、回転機械の動特性を演算している。

特公平３−６５８５７号公報特公平３−６５８５８号公報

しかしながら、特許文献１、２では、回転機械全体としての動特性を演算するに留まっており、回転機械のロータ系を構成する各機械要素（例えば回転軸、羽根車、軸受、シール部材など）の動特性を演算しているわけではない。したがって、特許文献１、２の技術を採用したところで、各機械要素の動特性は不明であり、ロータに自励振動が生じるか否かの判断、あるいはロータに自励振動が生じないようにするための設計の検討などが十分に行えないおそれがあった。

本発明は、上記課題に鑑みてなされたものであり、ターボ機械のロータ系を構成する各機械要素の動特性を演算するための技術を提供することを目的とする。

上記目的を達成するため、本発明におけるターボ機械の動特性演算方法は、回転軸に連結された羽根車を回転させることで流体を圧縮または圧送するターボ機械の動特性演算方法であって、前記ターボ機械のロータ系の振動方程式における特性行列を、事前に動特性が算出可能な機械要素に関する既知のパラメータと、事前に動特性が算出困難な機械要素に関する未知のパラメータとで表すことで、前記ロータ系を振動方程式によりモデル化するモデル化工程と、前記ロータ系に発生している振動を検出する振動検出工程と、前記振動検出工程で検出した振動の信号データに基づいて、前記ロータ系の固有値および振動モードを算出するモード解析工程と、前記モデル化工程で得た振動方程式と、前記モード解析工程で算出した固有値および振動モードとから、前記未知のパラメータを算出する動特性算出工程と、を備えることを特徴とする。

また、上記目的を達成するため、本発明におけるターボ機械の動特性演算装置は、回転軸に連結された羽根車を回転させることで流体を圧縮または圧送するターボ機械の動特性演算装置であって、前記ターボ機械のロータ系の振動方程式における特性行列を、事前に動特性が算出可能な機械要素に関する既知のパラメータと、事前に動特性が算出困難な機械要素に関する未知のパラメータとで表すことで、前記ロータ系を振動方程式によりモデル化するモデル化部と、前記ロータ系に発生している振動を検出する振動検出手段から出力された信号データに基づいて、前記ロータ系の固有値および振動モードを算出するモード解析部と、前記モデル化部で得た振動方程式と、前記モード解析部で算出した固有値および振動モードとから、前記未知のパラメータを算出する動特性算出部と、を備えることを特徴とする。

一般的に、ロータ系を構成する機械要素のうち、回転軸の剛性や質量は寸法や物性などから精度よく算出できるのに対し、羽根車、軸受、シール部材などの回転軸以外の機械要素の動特性を精度よく予測するのは難しい。そこで、本発明では、ターボ機械のロータ系の振動方程式における特性行列を、事前に動特性が算出可能な機械要素（例えば回転軸）に関する既知のパラメータと、事前に動特性が算出困難な機械要素（例えば羽根車、軸受、シール部材）に関する未知のパラメータとに分けて、ロータ系のモデル化を行っている。そして、このモデル化により得られた振動方程式と、ロータ系に発生している振動の信号データから算出されるロータ系の固有値および振動モードとから、上記未知のパラメータを算出している。すなわち、本発明によれば、上述のようにモデル化することで、動特性を予測しがたい機械要素に関する未知のパラメータも算出することができ、ターボ機械のロータ系を構成する各機械要素の動特性を演算可能となっている。

本発明にかかるターボ機械の動特性演算方法では、前記ロータ系を支持する土台系の動特性に関するパラメータを前記未知のパラメータに含めて前記モデル化を行うようにしてもよい。同様に、本発明にかかるターボ機械の動特性演算装置では、前記ロータ系を支持する土台系の動特性に関するパラメータを前記未知のパラメータに含めて前記モデル化を行うようにしてもよい。

ロータ系を支持する土台系の剛性がロータ系と比較して十分に大きければ、土台系の動特性の影響を無視しても特に支障はない。しかしながら、例えばターボ機械が大型でロータ系を支持する軸受などの剛性を大きくする必要がある場合には、相対的に土台系の剛性が小さくなる。このような場合には、上述のように、土台系の動特性に関するパラメータを未知のパラメータに含めて動特性の演算を行うことで、ロータ系を構成する各機械要素の動特性を精度よく演算することができる。

ここで、本発明にかかるターボ機械の動特性演算方法では、前記振動検出工程において、前記ターボ機械の外部に設けられた加振機で振動を加えてもよい。このように、加振機によりロータ系に加振できるように構成すれば、ロータ系の振動モードを励起させることができる。その結果、後のモード解析工程における固有値および振動モードの算出に適した振動データを測定することが可能となり、好適に各機械要素の動特性を演算することができる。

また、前記モード解析工程において、前記信号データを周波数分析して、前記ロータ系の固有値および振動モードを算出し、前記動特性算出工程において、最小二乗法により前記未知のパラメータを算出するとよい。このように、最小二乗法により未知のパラメータを算出することで、算出された値は常に最小二乗解であるため、安定して精度のよい値を算出することができる。その結果、各機械要素の動特性をより的確に把握することが可能となる。

本発明の第１実施形態にかかる動特性演算装置および測定対象であるターボ圧縮機を示す模式図である。図１のＡ−Ａ矢視図である。第１実施形態にかかる動特性演算方法を示すフローチャートである。第２実施形態にかかる動特性演算方法を示すフローチャートである。本発明の第３実施形態にかかる動特性演算装置および測定対象であるターボ圧縮機を示す模式図である。第３実施形態にかかる動特性演算方法を示すフローチャートである。第４実施形態にかかる動特性演算方法を示すフローチャートである。

以下、本発明にかかるターボ機械の動特性演算装置および動特性演算方法の実施形態について、図面を参照しつつ説明する。なお、以下の実施形態では、遠心式のターボ圧縮機を測定対象としているが、羽根車を回転させることで流体を圧縮または圧送するターボ機械であれば測定対象を他のものとしてもよい。例えば、軸流式のターボ圧縮機、ターボポンプ、ブロワ、タービンなどを測定対象とすることも可能である。

［第１実施形態］
図１、２を参照しつつ、本発明の第１実施形態にかかる動特性演算装置の構成について説明した後、図３を参照しつつ、本発明の第１実施形態にかかる動特性演算方法について説明する。

<ターボ圧縮機>
図１に示す遠心式のターボ圧縮機１は、動特性演算装置１００による動特性の測定対象となるものである。ターボ圧縮機１は、ハウジング１０の内部空間にロータ系２０が収容されて構成されている。ロータ系２０は、回転軸２１と、回転軸２１の両端部に設けられた一対の羽根車２２と、回転軸２１を回転自在に支持する２つの軸受２３と、羽根車２２と軸受２３との間に設けられた２つのシール部材２４とを有している。軸受２３としては、例えばすべり軸受を採用することができる。なお、ロータ系２０の具体的構成、すなわちどのような機械要素をいくつ、どこに配置するかは適宜変更が可能である。また、以下の説明においては、回転軸２１と一対の羽根車２２とで構成される回転体を、ロータ系２０とは区別して単にロータ２５と称する。

回転軸２１の中央部には、ギア（ピニオンギア）３０が設けられている。このギア３０は、ギア３０よりも大きな径の大径ギア（不図示）に噛み合わせられており、ロータ２５はギア３０を介して回転させられる。ギア３０と上記大径ギアとで増速機構が構成されており、例えば２００００〜３００００ｒｐｍの回転速度でロータ２５を回転させることができる。

ハウジング１０のうち、各羽根車２２の表面（回転軸２１とは反対側の面）と対向する部位には吸込口１１が形成されている。また、ハウジング１０の内壁によって画定される内部空間のうち、羽根車２２よりも径方向外側の部分は、羽根車２２が遠心方向に気体を吐き出すことができる吐出室１２として構成されている。このような構成により、ターボ圧縮機１は、図１の一点鎖線の矢印で示すように、羽根車２２により吸込口１１から気体を吸い込み、吸い込んだ気体を羽根車２２の回転により圧縮し、さらに圧縮した気体を吐出室１２に吐き出す。

<動特性演算装置>
次に、動特性演算装置１００について説明する。動特性演算装置１００は、ターボ圧縮機１のロータ系２０を構成する各機械要素（回転軸２１、羽根車２２、軸受２３、シール部材２４）の動特性を演算するための装置であり、モデル化部１０１、モード解析部１０２および動特性算出部１０３の各機能部を有して構成されている。これらの機能については、後述する。なお、ロータ系２０が他の機械要素を含む場合には、当該機械要素の動特性を演算することももちろん可能である。

動特性演算装置１００が有する上記各機能部は、マイクロコンピュータなどの論理演算ハードウェアを中核とし、当該ハードウェア上で実行されるプログラムなどのソフトウェアとの協働により実現される。すなわち、各機能部は、機能としての分担を示すものであり、必ずしも物理的に独立して構成される必要はない。なお、プログラムは不図示のメモリに記憶されており、動特性演算装置１００が必要に応じてメモリからプログラムを読み出すことで、各機械要素の動特性の演算が行われる。

<入力手段および表示手段>
動特性演算装置１００には、入力手段１１１および表示手段１１２が接続されている。入力手段１１１は、動特性演算装置１００に対してデータの入力を行うための手段であり、キーボードなどにより構成されている。作業者は、入力手段１１１を介して、例えば後述する既知のパラメータを動特性演算装置１００に入力することができる。

また、表示手段１１２は、動特性演算装置１００から出力されたデータを作業者が認識できるように表示するための手段であり、モニタなどにより構成されている。例えば、後述する未知のパラメータを動特性演算装置１００により算出した結果を、表示手段１１２に表示させることができる。

<振動センサー>
ターボ圧縮機１の回転軸２１の周面近傍には、回転軸２１の振動を検出するための振動センサー（振動検出手段）３１が複数配設されている。振動センサー３１は、回転軸２１の軸方向に関しては、図１に示すように４箇所に設けられており、回転軸２１の周方向に関しては、図２に示すように２箇所に設けられている。なお、図２は図１のＡ−Ａ矢視図である。

より具体的には、軸方向においては、羽根車２２とシール部材２４との間の位置、および軸受２３とシール部材２４との間の位置にそれぞれ振動センサー３１が設けられている。さらに、軸方向の各位置においては、図２に示すように、回転軸２１の上部の周面に沿うように、例えば９０°の所定間隔を空けて、鉛直面に対して対称に２つずつ配置されている。ただし、このような振動センサー３１の配置形態は一例にすぎず、その個数や位置は適宜変更が可能である。

振動センサー３１は、例えば回転軸２１の変位を検出することで、回転軸２１の振動に関する信号データを生成し、これを動特性演算装置１００に出力する。振動センサー３１が出力する回転軸２１の振動に関する信号データは、後述する機械要素の動特性の演算において、ロータ系２０の振動に関するデータとして用いられる。

<加振機>
本実施形態では、ターボ圧縮機１に振動を加えるための加振機３２を設けている。加振機３２は、例えば電磁式のものであり、ターボ圧縮機１のハウジング１０上に載置されている。加振機３２を作動させると、加振機３２で発生した振動は、ハウジング１０、軸受２３を経由して、間接的に回転軸２１に伝わる。なお、加振機３２で加えることのできる振動の条件は、適宜変更することが可能である。

<動特性の演算>
動特性演算装置１００によるロータ系２０の各機械要素の動特性の演算方法を、図３を参照しつつ説明する。ここでは、例として、事前に動特性が算出可能な機械要素として回転軸２１を取り上げ、事前に動特性が算出困難な機械要素として軸受２３を取り上げる。そして、軸受２３の動特性として、軸受２３の剛性行列Ｋ_brおよび減衰行列Ｃ_brを求める場合について説明する。なお、以下の各数式では、行列を太字の大文字で示し、ベクトルを太字の小文字で示しているが、明細書の文中においては太字が使用できないため、通常の字体で記載している。

（モデル化工程）
ここで、λ_iをロータ系２０のｉ次の固有値、φ_iをｉ次の固有値に対応する振動モードとすると、軸受２３を含むロータ系２０の振動方程式は、ロータ系２０の特性行列（質量行列Ｍ、減衰行列Ｃおよび剛性行列Ｋ）を用いて式（１）により表すことができる。

そして、本実施形態にかかる動特性の演算方法においては、まず最初にロータ系２０のモデル化を行う（ステップＳ２０１）。ここで言うモデル化とは、ロータ系２０の振動方程式（１）における特性行列Ｍ、Ｃ、Ｋを、回転軸２１に関する既知のパラメータと、軸受２３に関する未知のパラメータとで表すことを指す。

具体的には、設計図面などから得られる寸法や物性の情報に基づいて、事前に算出可能な既知のパラメータ、すなわち回転軸２１の質量行列Ｍ_rt、減衰行列Ｃ_rtおよび剛性行列Ｋ_rtについては、あらかじめ作業者がこれらの値を算出し、その値を入力手段１１１を介して動特性演算装置１００に入力する。一方、事前に算出することが困難な未知のパラメータ、すなわち軸受２３の剛性行列Ｋ_brおよび減衰行列Ｃ_brについては、これらを未知数としておく。

既知のパラメータと未知のパラメータとを用いて式（１）を書き改めると、式（２）のようになる。そして、式（２）の未知のパラメータを含む項を左辺に移項し、既知のパラメータを右辺に移項して式を整理すると式（３）が得られる。

このように、ロータ系２０の振動方程式における特性行列を、既知のパラメータと未知のパラメータとで表すモデル化の処理は、動特性演算装置１００のモデル化部１０１にて実行される。

（振動検出工程）
次に、ターボ圧縮機１の運転を開始し（ステップＳ２０２）、続いて加振機３２を作動しターボ圧縮機１を加振する（ステップＳ２０３）。このとき、回転軸２１は加振機３２により間接的に加振される。そして、ロータ系２０に発生している振動を、複数の振動センサー３１により検出する（ステップＳ２０４）。なお、加振機３２による加振形態としては種々の形態を採用し得るが、ここではスイープ加振を行うものとする。

（モード解析工程）
動特性演算装置１００のモード解析部１０２は、振動センサー３１から出力された信号データに基づいて、ロータ系２０の固有振動数、減衰比および振動モードを算出するモード解析を実行する（ステップＳ２０５）。具体的には、スイープ加振が行われたロータ系２０の周波数応答関数を測定し、例えばＭＤＯＦ法を用いることで上述の各値を算出することができる。こうして、式（３）におけるロータ系２０の固有値λ_i（固有振動数および減衰比から算出）および振動モードφ_iが求められる。

（動特性算出工程）
モード解析工程において、ロータ系２０の固有値λ_iおよび振動モードφ_iが求められると、式（３）の右辺は既知の値となる。そこで、式（３）の右辺を既知ベクトルｂとし、左辺を各成分に展開すると、式（４）が得られる。

ここで、行列成分Ｃ、Ｋおよび振動モードφの添え字のうち、「１」は１つ目の軸受２３に関するものであること、「２」は２つ目の軸受２３に関するものであることを示す。また、行列成分Ｃ、Ｋの添え字のうち、「ｘｘ」は水平方向に変位した場合に生じる水平方向の減衰力または弾性力に関するものであること、「ｘｙ」は水平方向に変位した場合に生じる垂直方向の減衰力または弾性力に関するものであること、「ｙｘ」は垂直方向に変位した場合に生じる水平方向の減衰力または弾性力に関するものであること、「ｙｙ」は垂直方向に変位した場合に生じる垂直方向の減衰力または弾性力に関するものであることを示す。また、振動モードφの添え字のうち、「ｘ」は水平方向に関する振動モード、「ｙ」は垂直方向に関する振動モードであることを示す。なお、軸受２３に関する振動モードφは、振動センサー３１から出力された信号データに基づいて求められる。

式（４）の左辺を既知の係数行列と未知のベクトルとに整理すると、式（５）が得られる。

ただし、

である。

さらに、式（５）の係数行列をＡとし、実部と虚部とに分離すると、式（７）が得られる。

式（７）は各モードに対してそれぞれ成り立つので、式（８）に示すように１次〜ｍ次モードまで使用する。

式（８）の係数行列の行数が列数より多い場合、すなわち未知のパラメータの数より式の数が多ければ、最小二乗法で解くことができる。このときの解を式（９）に示す。ただし、添え字の「Ｔ」は転置行列を、「−１」は逆行列を表す。

このように、モデル化工程で得た振動方程式（３）と、モード解析工程で算出した固有値λ_iおよび振動モードφ_iとから、未知のパラメータ、すなわち軸受２３の剛性行列Ｋ_brおよび減衰行列Ｃ_brを算出する処理は、動特性演算装置１００の動特性算出部１０３にて実行される（ステップＳ２０６）。

最後に、動特性演算装置１００は、測定を終了するかどうかを判断し（ステップＳ２０７）、測定を終了する場合には動作を終了し、測定を終了しない場合にはステップＳ２０３に戻って、以降のステップを繰り返し実行する。

動特性算出部１０３により算出された軸受２３の剛性行列Ｋ_brおよび減衰行列Ｃ_brに関する情報は、測定終了時、または測定中の適宜のタイミングで、表示手段１１２に表示される。そして、表示手段１１２に表示された結果を参照することで、作業者はロータ２５に自励振動が生じるか否かを判断したり、あるいはロータ２５に自励振動が生じないようにロータ系２０の各機械要素の設計の検討を適切に行うことが可能となる。

以上のように、本実施形態では、ターボ圧縮機１のロータ系２０の振動方程式（１）における特性行列Ｍ、Ｃ、Ｋを、事前に動特性が算出可能な機械要素である回転軸２１に関する既知のパラメータ（回転軸２１の質量行列Ｍ_rt、減衰行列Ｃ_rtおよび剛性行列Ｋ_rt）と、事前に動特性が算出困難な機械要素である軸受２３に関する未知のパラメータ（軸受２３の剛性行列Ｋ_brおよび減衰行列Ｃ_br）とに分けて、ロータ系２０のモデル化を行っている。そして、このモデル化により得られた振動方程式（３）と、ロータ系２０に発生している振動の信号データから算出されるロータ系２０の固有値λ_iおよび振動モードφ_iとから、上記未知のパラメータを算出している。すなわち、本実施形態によれば、上述のようにモデル化することで、動特性を予測しがたい機械要素（軸受２１など）に関するパラメータも算出することができ、ターボ圧縮機１のロータ系２０を構成する各機械要素の動特性を演算可能となっている。

また、本実施形態では、振動検出工程において、ターボ圧縮機１の外部に設けられた加振機３２で振動を加えている。このように、加振機３２によりロータ系２０に加振できるように構成すれば、ロータ系２０の振動モードを励起させることができる。その結果、後のモード解析工程における固有値および振動モードの算出に適した振動データを測定することが可能となり、好適に各機械要素の動特性を演算することができる。

また、本実施形態では、モード解析工程において、振動センサー３１の出力した信号データを周波数分析して（周波数応答関数を測定して）、ロータ系２０の固有値および振動モードを算出し、動特性算出工程において、最小二乗法により未知のパラメータを算出している。このように、最小二乗法により未知のパラメータを算出することで、算出された値は常に最小二乗解であるため、安定して精度のよい値を算出することができる。その結果、各機械要素の動特性をより的確に把握することが可能となる。

［第２実施形態］
次に、図４を参照しつつ、本発明の第２実施形態にかかる動特性演算方法について説明する。本実施形態における動特性演算方法が、第１実施形態の動特性演算方法と異なる点は、ターボ圧縮機１を加振機３２で加振する工程が省略されている点である。その他の点については第１実施形態と同様であり、図４に示すステップＳ３０１〜Ｓ３０６は、図３に示すステップＳ２０１、Ｓ２０２、Ｓ２０４〜Ｓ２０７と同様の工程であるので、これらの工程についての説明は省略する。

本実施形態のように、加振機３２による加振を行わなくても、ロータ２５にはターボ圧縮機１の運転に伴う振動が発生するため、動特性の測定は可能である。そこで、加振機３２により加振する工程を省くことで、加振機３２を設置するコストや手間が省けるため、より簡単にロータ系２０の各機械要素の動特性の測定を行うことができる。

［第３実施形態］
次に、図５および図６を参照しつつ、本発明の第３実施形態にかかる動特性演算方法について説明する。本実施形態における動特性演算方法が、第１実施形態の動特性演算方法と異なる点は、ロータ系２０を支持する土台系４０に関する動特性も考慮して、ロータ系２０に関する運動方程式を解いていく点にある。なお、第１実施形態と同じ構成要素に対しては同じ符号を付し、説明を省略する。以下、第１実施形態と異なる点を中心に説明する。

図５に示すように、ターボ圧縮機１は、ロータ系２０を収容するハウジング１０が、支持部材４１を介して固定面Ｇに対して固定支持されている。以下、ハウジング１０と支持部材４１とを合わせて土台系４０と称する。ここで、土台系４０の剛性がロータ系２０に対して十分に大きい場合は、第１実施形態でそうしたように、土台系４０の動特性の影響を無視してロータ系２０の動特性を演算しても特に問題はない。

しかしながら、例えばターボ圧縮機１が大型となると、それに伴ってロータ系２０を支持する軸受２３の剛性を大きくする必要があり、その結果、相対的に土台系４０の剛性が小さくなることがある。そうすると、土台系４０の動特性の影響を無視することができず、ロータ系２０の動特性を精度よく演算するためには、土台系４０の動特性を考慮する必要が生じる。本実施形態は、このようなケースを想定して、土台系４０の動特性も考慮して、ロータ系２０の動特性を演算するものである。

図５に示すターボ圧縮機１では、回転軸２１の中央部に設けられたギア３０に対して、ギア３０よりも大きな径の大径ギア３３が噛み合わせられて、増速機構が構成されている。そして、大径ギア３３を収容する空間を確保するため、ハウジング１０の中央部が一部外側に突出して、突出部１０ａが形成されている。上記支持部材４１は、この突出部１０ａの左右両側に設けられている。すなわち、ハウジング１０および支持部材４１からなる土台系４０は、突出部１０ａを中心として概ね左右対称に構成されている。なお、加振機３２は突出部１０ａ上に配設される。

後で説明する動特性の演算においては、土台系４０を、突出部１０ａよりも左側のハウジング１０の部分と左側の支持部材４１からなる系と、突出部１０ａよりも右側のハウジング１０の部分と右側の支持部材４１からなる系との２つに分けてモデル化している。ただし、土台系４０のモデル化方法はこれに限定されず、実際の構成に適合するように適宜変更が可能である。

本実施形態においては、回転軸２１の振動を検出するための振動センサー３１に加えて、土台系４０の振動を検出するための振動センサー３４が複数配設されている。具体的には、振動センサー３４は、ハウジング１０の外側、且つ軸受２３と対向する位置に設けられ、ハウジング１０の変位を検出することで、土台系４０の振動に関する信号データを生成して、動特性演算装置１００に出力する。このように、軸受２３の近傍に振動センサー３４を設けることで、軸受２３を介して土台系４０からロータ系２０へ伝播される振動を精度よく検出することができる。ただし、振動センサー３４の位置や個数は適宜変更が可能である。

<動特性の演算>
次に、図６を参照しつつ、動特性演算装置１００によるロータ系２０の各機械要素の動特性の演算方法について説明する。ここでは、例として、事前に動特性が算出可能な機械要素として回転軸２１を取り上げ、事前に動特性が算出困難な機械要素として軸受２３を取り上げ、さらに事前に算出が困難な土台系４０の動特性についても考慮する。そして、軸受２３の動特性として、軸受２３の剛性行列Ｋ_brおよび減衰行列Ｃ_brを求めるとともに、土台系４０の動特性として、土台系４０の質量行列Ｍ_bsおよび剛性行列Ｋ_bsを求める場合について説明する。なお、以下の各数式では、行列を太字の大文字で示し、ベクトルを太字の小文字で示しているが、明細書の文中においては太字が使用できないため、通常の字体で記載している。

（モデル化工程）
ここで、λ_iをロータ系２０のｉ次の固有値、φ_iをｉ次の固有値に対応する振動モードとすると、軸受２３と土台系４０を含むロータ系２０の振動方程式は、ロータ系２０の特性行列（質量行列Ｍ、減衰行列Ｃおよび剛性行列Ｋ）を用いて式（１０）により表すことができる。

そして、本実施形態にかかる動特性の演算方法においては、まず最初に土台系４０を含めたロータ系２０のモデル化を行う（ステップＳ４０１）。ここで言うモデル化とは、ロータ系２０の振動方程式（１０）における特性行列Ｍ、Ｃ、Ｋを、回転軸２１に関する既知のパラメータと、軸受２３および土台系４０に関する未知のパラメータとで表すことを指す。

具体的には、設計図面などから得られる寸法や物性の情報に基づいて、事前に算出可能な既知のパラメータ、すなわち回転軸２１の質量行列Ｍ_rt、減衰行列Ｃ_rtおよび剛性行列Ｋ_rtについては、あらかじめ作業者がこれらの値を算出し、その値を入力手段１１１を介して動特性演算装置１００に入力する。一方、事前に算出することが困難な未知のパラメータ、すなわち軸受２３の剛性行列Ｋ_brおよび減衰行列Ｃ_br、ならびに土台系４０の質量行列Ｍ_bsおよび剛性行列Ｋ_bsについては、これらを未知数としておく。なお、ここでは、土台系４０の減衰は、ロータ系２０と比較して十分に小さいものとして考慮していない。

既知のパラメータと未知のパラメータとを用いて式（１０）を書き改めると、式（１１）のようになる。そして、式（１１）の未知のパラメータを含む項を左辺に移項し、既知のパラメータを右辺に移項して式を整理すると式（１２）が得られる。

このように、土台系４０の動特性を考慮したロータ系２０の振動方程式における特性行列を、既知のパラメータと未知のパラメータとで表すモデル化の処理は、動特性演算装置１００のモデル化部１０１にて実行される。

（振動検出工程）
次に、ターボ圧縮機１の運転を開始し（ステップＳ４０２）、続いて加振機３２を作動しターボ圧縮機１を加振する（ステップＳ４０３）。このとき、回転軸２１は加振機３２により間接的に加振される。そして、ロータ系２０に発生している振動を複数の振動センサー３１により検出するとともに、土台系４０に発生している振動を複数の振動センサー３４により検出する（ステップＳ４０４）。なお、加振機３２による加振形態としては種々の形態を採用し得るが、ここではスイープ加振を行うものとする。

（モード解析工程）
動特性演算装置１００のモード解析部１０２は、振動センサー３１から出力された信号データに基づいて、ロータ系２０の固有振動数、減衰比および振動モードを算出するモード解析を実行する（ステップＳ４０５）。具体的には、スイープ加振が行われたロータ系２０の周波数応答関数を測定し、例えばＭＤＯＦ法を用いることで上述の各値を算出することができる。こうして、式（１２）におけるロータ系２０の固有値λ_i（固有振動数および減衰比から算出）および振動モードφ_iが求められる。

（動特性算出工程）
モード解析工程において、ロータ系２０の固有値λ_iおよび振動モードφ_iが求められると、式（１２）の右辺は既知の値となる。そこで、式（１２）の右辺を既知ベクトルｂとし、左辺を各成分に展開すると、式（１３）が得られる。

ここで、行列成分Ｍ、Ｃ、Ｋおよび振動モードφの添え字のうち、「ｂｒ１」は図５の左側の軸受２３に関するものであること、「ｂｒ２」は図５の右側の軸受２３に関するものであること、「ｂｓ１」は図５の左側半分の土台系４０に関するものであること、「ｂｓ２」は図５の右側半分の土台系４０に関するものであることを示す。また、行列成分Ｍ、Ｃ、Ｋの添え字のうち、「ｘｘ」は水平方向に変位した場合に生じる水平方向の慣性力、減衰力または弾性力に関するものであること、「ｘｙ」は水平方向に変位した場合に生じる垂直方向の慣性力、減衰力または弾性力に関するものであること、「ｙｘ」は垂直方向に変位した場合に生じる水平方向の慣性力、減衰力または弾性力に関するものであること、「ｙｙ」は垂直方向に変位した場合に生じる垂直方向の慣性力、減衰力または弾性力に関するものであることを示す。また、振動モードφの添え字のうち、「ｘ」は水平方向に関する振動モード、「ｙ」は垂直方向に関する振動モードであることを示す。なお、軸受２３に関する振動モードφは、振動センサー３１から出力された信号データに基づいて求められ、土台系４０に関する振動モードφは、振動センサー３４から出力された信号データに基づいて求められる。

式（１３）の左辺を既知の係数行列と未知のベクトルとに整理すると、式（１４）が得られる。

ここで、Φ_i,M、Φ_i,C、Φ_i,Kは、式（１３）の左辺から未知数であるＭ、Ｃ、Ｋをベクトルとしてまとめ、かつ、スカラー量であるλをくくり出して残った成分を示す。また、

である。なお、式（１５）において軸受２３と土台系４０とは同じ式であるので、添え字の「ｂｒ」と「ｂｓ」は省略している。

さらに、式（１４）の係数行列をＡとし、実部と虚部とに分離すると、式（１６）が得られる。

式（１６）は各モードに対してそれぞれ成り立つので、式（１７）に示すように１次〜ｍ次モードまで使用する。

式（１７）の係数行列の行数が列数より多い場合、すなわち未知のパラメータの数より式の数が多ければ、最小二乗法で解くことができる。このときの解を式（１８）に示す。ただし、添え字の「Ｔ」は転置行列を、「−１」は逆行列を表す。

このように、モデル化工程で得た振動方程式（１２）と、モード解析工程で算出した固有値λ_iおよび振動モードφ_iとから、未知のパラメータ、すなわち軸受２３の剛性行列Ｋ_brおよび減衰行列Ｃ_br、ならびに土台系４０の質量行列Ｍ_bsおよび剛性行列Ｋ_bsを算出する処理は、動特性演算装置１００の動特性算出部１０３にて実行される（ステップＳ４０６）。

最後に、動特性演算装置１００は、測定を終了するかどうかを判断し（ステップＳ４０７）、測定を終了する場合には動作を終了し、測定を終了しない場合にはステップＳ４０３に戻って、以降のステップを繰り返し実行する。

動特性算出部１０３により算出された軸受２３の剛性行列Ｋ_brおよび減衰行列Ｃ_br、ならびに土台系４０の質量行列Ｍ_bsおよび剛性行列Ｋ_bsに関する情報は、測定終了時、または測定中の適宜のタイミングで、表示手段１１２に表示される。そして、表示手段１１２に表示された結果を参照することで、作業者はロータ２５に自励振動が生じるか否かを判断したり、あるいはロータ２５に自励振動が生じないようにロータ系２０の各機械要素や土台系４０の設計の検討を適切に行うことが可能となる。

以上のように、本実施形態では、ターボ圧縮機１のロータ系２０の振動方程式（１０）における特性行列Ｍ、Ｃ、Ｋを、事前に動特性が算出可能な機械要素である回転軸２１に関する既知のパラメータ（回転軸２１の質量行列Ｍ_rt、減衰行列Ｃ_rtおよび剛性行列Ｋ_rt）と、事前に動特性が算出困難な機械要素である軸受２３および土台系４０に関する未知のパラメータ（軸受２３の剛性行列Ｋ_brおよび減衰行列Ｃ_br、ならびに土台系４０の質量行列Ｍ_bsおよび剛性行列Ｋ_bs）とに分けて、土台系４０を含めたロータ系２０のモデル化を行っている。そして、このモデル化により得られた振動方程式（１２）と、ロータ系２０に発生している振動の信号データから算出されるロータ系２０の固有値λ_iおよび振動モードφ_iとから、上記未知のパラメータを算出している。すなわち、本実施形態によれば、上述のようにモデル化することで、動特性を予測しがたい機械要素（軸受２１など）に関するパラメータも算出することができ、ターボ圧縮機１のロータ系２０を構成する各機械要素の動特性を演算可能となっている。しかも、土台系４０の動特性に関するパラメータを未知のパラメータに含めて動特性の演算を行うことで、土台系４０の動特性を無視できない場合においても、ロータ系２０を構成する各機械要素の動特性を精度よく演算することができる。

［第４実施形態］
次に、図７を参照しつつ、本発明の第４実施形態にかかる動特性演算方法について説明する。本実施形態における動特性演算方法が、第３実施形態の動特性演算方法と異なる点は、ターボ圧縮機１を加振機３２で加振する工程が省略されている点である。その他の点については第３実施形態と同様であり、図７に示すステップＳ５０１〜Ｓ５０６は、図６に示すステップＳ４０１、Ｓ４０２、Ｓ４０４〜Ｓ４０７と同様の工程であるので、これらの工程についての説明は省略する。

［その他］
本発明は上記実施形態に限定されるものではなく、その趣旨を逸脱しない限りにおいて上記実施形態の要素を適宜組み合わせまたは種々の変更を加えることが可能である。

上記実施形態では、事前に動特性が算出可能な機械要素として回転軸２１を取り上げ、事前に動特性が算出困難な機械要素として軸受２３を取り上げたが、どの機械要素の動特性を既知のパラメータあるいは未知のパラメータとして扱うかは、適宜変更が可能である。例えば、羽根車２２やシール部材２４などの回転軸２１以外の機械要素に関する動特性のパラメータを未知数とすることで、羽根車２２やシール部材２４に関する動特性を演算することも可能である。

また、上記実施形態では、ロータ系２０のモデル化（ステップＳ２０１、Ｓ３０１）を、ターボ圧縮機１の運転を開始（ステップＳ２０２、Ｓ３０２）する前に実行するものとしたが、ターボ圧縮機１の運転を開始した後にモデル化を実行するようにしてもよい。

１ターボ圧縮機（ターボ機械）
１０ハウジング
２０ロータ系
２１回転軸
２２羽根車
２３軸受
２４シール部材
３１振動センサー（振動検出手段）
３２加振機
３４振動センサー
４０土台系
４１支持部材
１００動特性演算装置
１０１モデル化部
１０２モード解析部
１０３動特性算出部

Claims

回転軸に連結された羽根車を回転させることで流体を圧縮または圧送するターボ機械の動特性演算方法であって、
前記ターボ機械のロータ系の振動方程式における特性行列を、事前に動特性が算出可能な機械要素に関する既知のパラメータと、事前に動特性が算出困難な機械要素に関する未知のパラメータとで表すことで、前記ロータ系を振動方程式によりモデル化するモデル化工程と、
前記ロータ系に発生している振動を検出する振動検出工程と、
前記振動検出工程で検出した振動の信号データに基づいて、前記ロータ系の固有値および振動モードを算出するモード解析工程と、
前記モデル化工程で得た振動方程式と、前記モード解析工程で算出した固有値および振動モードとから、前記未知のパラメータを算出する動特性算出工程と、
を備えることを特徴とするターボ機械の動特性演算方法。
前記ロータ系を支持する土台系の動特性に関するパラメータを前記未知のパラメータに含めて前記モデル化を行う請求項１に記載のターボ機械の動特性演算方法。
前記振動検出工程において、前記ターボ機械の外部に設けられた加振機で振動を加える請求項１または２に記載のターボ機械の動特性演算方法。
前記モード解析工程において、前記信号データを周波数分析して、前記ロータ系の固有値および振動モードを算出し、
前記動特性算出工程において、最小二乗法により前記未知のパラメータを算出する請求項１ないし３のいずれか１項に記載のターボ機械の動特性演算方法。
回転軸に連結された羽根車を回転させることで流体を圧縮または圧送するターボ機械の動特性演算装置であって、
前記ターボ機械のロータ系の振動方程式における特性行列を、事前に動特性が算出可能な機械要素に関する既知のパラメータと、事前に動特性が算出困難な機械要素に関する未知のパラメータとで表すことで、前記ロータ系を振動方程式によりモデル化するモデル化部と、
前記ロータ系に発生している振動を検出する振動検出手段から出力された信号データに基づいて、前記ロータ系の固有値および振動モードを算出するモード解析部と、
前記モデル化部で得た振動方程式と、前記モード解析部で算出した固有値および振動モードとから、前記未知のパラメータを算出する動特性算出部と、
を備えることを特徴とするターボ機械の動特性演算装置。
前記ロータ系を支持する土台系の動特性に関するパラメータを前記未知のパラメータに含めて前記モデル化を行う請求項５に記載のターボ機械の動特性演算装置。