レイリープリセットの式を解いてみる

 

はじめに

 

今日はレイリープリセットの式の数値解析を行っていきたいと思います。そもそも、この式は何かと申し上げますと、一言で気泡振動の運動方程式です。(内部の不揮発性気体を理想気体として扱うetc...等の条件は何かと多いですが)matlabのode45という汎用性のあるソルバで解析を行います。ode45は小田急沿線のアイドルグループではなく、常微分方程式の近似解を求めるためにmathworks社が用意してくれている関数です。ode45はルンゲクッタ法を採用しているらしいですが、詳しい中身はあんま知らないです。オイラー法しか知りません。きっとオイラー法より精度がいいんでしょうか?とりあえず、用意して下さったので少しばかり拝借させて頂きます。

 

f:id:yama3394145:20170112020353p:plain

レイリー・プレセット方程式 - Wikipedia

 

やってみる

 

まず、関数ファイルRP.mを作成します。

f:id:yama3394145:20170112012419p:plain

また、入力する周波数のパラメータを弄りやすいように、pulse.mも作成しとく

f:id:yama3394145:20170112013205p:plain

これ解きたい微分方程式の定義が完了したので、これをode45に素直にぶち込む。

はじめは、R0=1[µm],として音圧の周波数f[Hz]を変化させて試してみる。図は時間変化に対する無次元数(Rを初期気泡径R0で割った)で表した。

f:id:yama3394145:20170112014205p:plain

こちらはf=500kHzのもの。1.5µs以降は、基本同じ周期で振動を繰り返す。ほんの初期の0.2µs付近では急激に落ち込んでいる。気泡体積も初期の1/8になっている。断熱圧縮を仮定して温度や圧力を計算するのかな?そこら辺は追々考えを纏めます。

f:id:yama3394145:20170112015733p:plain

続いて、f=1MHz同様の気泡径で検証。とりあえず微分方程式が比較的簡単に解けることが分かったので良し。後は条件を変えてやる計算をひたすら延々と繰り返す。一番の問題はどう評価するか?これはもう少し気泡力学の知識を増やす必要がある。