いきなりOpenFOAM (65)

インパルス応答解析の自動化(その2)

csvファイルの読み取り

 今回は、いきなりOpenFOAM第63回で計算したタンク内のスロッシングのインパルス応答波形から振幅倍率の周波数応答特性を求めるコードを作成してみます。モデルの詳細や解析条件は、いきなりOpenFOAM第63回を参照してください。液面の時間変動はcsvファイルとして出力され、その内の安定した部分を取り出したheight.csvファイルから振幅倍率を求めます。
 csvファイルの読取りはnumpy.loadtxt(ファイル名)で行います。読取った配列の1列目を時間t、2列目を液面高さyとして、yは平均値からの差分に変換します。以上のコードが図1です。念のため、読み取った結果を図として表示させたものが図2です。

図1 csvファイルの読み取り
図2 読み取ったcsvファイルの波形
FFT処理と固有角周波数および減衰比の算出

 次に、FFT処理を行うには、時間刻みが一定でなければなりません。そこで、scipyのinterpolateモジュールの一次元補間interp1dモジュールを用いて、時間刻みが一定のデータを作成します。具体的には図3に示すように、関数をinterpolate.interp1d(参照xデータ、参照yデータ、kind=補間の種類)として設定して、この関数にarange(40,60,0.01)で設定した0.01秒刻みのデータを与えて補間された液面の時間変動を得ます。

図3 均一時間刻みのデータに補間

データの時間刻みが一定となったので、いきなりOpenFOAM第60回で説明した手法を用いて、図4に示すコードで固有角周波数と減衰比を求めます。

図4 固有角周波数と減衰比を計算
プログラムコード全体と動作結果

 以上をまとめると、図5示すコードとなり、これを動作させると、図6に示す振幅倍率の周波数特性が得られます。また、図6の結果は、numpy.savetextでMd.csvというcsvファイルに出力されます。これを表計算ソフトで読取り、振幅倍率の周波数特性を表示した結果が図7の橙色の曲線です。図中、青色の曲線は、いきなりOpenFOAM第62回で求めた隣り合った極大値同士の値から求めた振幅倍率です。図7から、一組のデータから求めた振幅倍率と全体から求めた振幅倍率とはほぼ同じ結果となることがわかります。

図5 読み取ったcsvファイルから振幅倍率の周波数特性を計算
図6 振幅倍率の周波数特性
図7 一組のデータから求めた振幅倍率と全体から求めた振幅倍率

 今回は、実際のスロッシング解析の結果から振幅倍率を求めました。実際のタンクでは、適当な仕切り板を設けることで、液体の揺動を抑えています。
 次回は、タンク内に仕切り板を設けた場合のスロッシング解析を行い、そのインパルス応答解析を行ってみます。

 このページでは、各アプリケーションの操作説明は省略しています。FreeCADの具体的な操作については、いきなりOpenFOAM第5回および第7回、OpenFOAMでの計算実行は第8回、ParaViewの操作については第3回第4回および第8回を参考にしてみてください。

おことわり
 本コンテンツの動作や表示はお使いのバージョンにより異なる場合があります。
 本コンテンツの動作ならびに設定項目等に関する個別の情報提供およびサポートはできかねますので、あらかじめご了承ください。
 本コンテンツは動作および結果の保証をするものではありません。ご利用に際してはご自身の判断でお使いいただきますよう、お願いいたします。