038_EvactChamber01.mw

チャンバの排気

イントロダクション

このワークシートでは、チューブから真空ポンプで接続された容積 V のチャンバに排気される動作を考えます。

チャンバから排気される際の圧力 P について計算を行います。

[Inserted Image]

> restart;

問題の記述

各パラメータは、容積 V は Float(5, -1)*m^2 l の長さが m 、直径 d が Float(5, -2)*m 、気温 0 ℃でのガスの粘度は Float(172, -1)*10^(-6)*Pa 、ポンプの排気の速度 S は 100*m^3/h 、チャンバの圧力は p で、真空ポンプの圧力は p0 となり初期圧力は Float(1013, -1)*kPa となります。

> val:= V=0.5, l=1, d=0.05, eta=17.2*10^(-6), S=100/3600:

理論

チャンバの圧力 p は以下の微分方程式を満たします。ここで、S_eff はチャンバのポンプの排気の速度です:

> deq := V*diff(p(t),t) = -p(t)*S_eff;

deq := V*(diff(p(t), t)) = -p(t)*S_eff

始点に付近の流れは連続関数になります。従って:

> eq1 := p(t)*S_eff = p0(t)*S;

eq1 := p(t)*S_eff = p0(t)*S

p0 がポンプの圧力であり、S はポンプの排気率になります。


チューブの薄板から成る流れのコンダクタンスを
以下で与えます:

> eqL := L = (Pi*(d/2)^4*(p(t)+p0(t)))/(8*eta*l*2);

eqL := L = 1/256*Pi*d^4*(p(t)+p0(t))/(eta*l)

また、違った圧力をかけた際の流れの割合の比率を表現します:

> eq2 := L = (p(t)*S_eff)/(p(t)-p0(t));

eq2 := L = p(t)*S_eff/(p(t)-p0(t))

解決法

これらの方程式を一つの微分方程式にまとめ定義します。

まず、eq1 を S_eff について解きます。

> tmp1:=solve(eq1,{S_eff});

tmp1 := {S_eff = p0(t)*S/p(t)}

eqL と eq2 を等式とし S_eff を代入し p0 について方程式の解を求めます。

> subs(eqL,tmp1,eq2);

1/256*Pi*d^4*(p(t)+p0(t))/(eta*l) = p0(t)*S/(p(t)-p0(t))

> tmp2:=solve(%, p0(t));

tmp2 := -(128*S*eta*l-(16384*S^2*eta^2*l^2+Pi^2*d^8*p(t)^2)^(1/2))/(Pi*d^4), -(128*S*eta*l+(16384*S^2*eta^2*l^2+Pi^2*d^8*p(t)^2)^(1/2))/(Pi*d^4)

2 つの解のうちの 1 つの正の p0 の値が物理的に正確な値となります。

これらの式の中にパラメータを代入し、解が正になるか確認します。

> evalf(subs(val, p(t) = 101.3*10^3, tmp2[1]));

101296.8854

> evalf(subs(val, p(t) = 101.3*10^3, tmp2[2]));

-101303.1147

従って:

> tmp3 := tmp2[1];

tmp3 := -(128*S*eta*l-(16384*S^2*eta^2*l^2+Pi^2*d^8*p(t)^2)^(1/2))/(Pi*d^4)

もしくはこの操作を自動的に出来るようにします:

> tmp3 := p0(t) = select(x -> evalb(0<evalf(subs(val, p(t)=101.3*10^3,x))),[tmp2])[1];

tmp3 := p0(t) = -(128*S*eta*l-(16384*S^2*eta^2*l^2+Pi^2*d^8*p(t)^2)^(1/2))/(Pi*d^4)

微分方程式に S_eff と p0 をそれぞれ代入し非線型微分方程式を得ます:

> DEQ:=subs(tmp1,tmp3,deq);

DEQ := V*(diff(p(t), t)) = (128*S*eta*l-(16384*S^2*eta^2*l^2+Pi^2*d^8*p(t)^2)^(1/2))*S/(Pi*d^4)

微分方程式にパラメータを入力し、数値的に解を求めます。この時、初期条件として、時間 t = 0、チャンバの標準空気圧を 101.3kPa として与えます。

> sol := dsolve({subs(val, DEQ), p(0)=101.3*10^3}, p(t), type=numeric);

sol := proc (x_rkf45) local res, data, vars, solnproc, outpoint, ndsol, i; option `Copyright (c) 2000 by Waterloo Maple Inc. All rights reserved.`; _EnvDSNumericSaveDigits := Digits; Digits := 14; if ...

このプロシージャ(微分方程式の解のプロシージャ)を使用し、ポンプで使用し始めた100秒後のチャンバの圧力を求めることが出来ます。

> sol(100);

[t = 100., p(t) = 394.713613604430122]

上記の100秒後の圧力は、sol の2つ目の結果の右辺になります。

また、下記の rhs(sol(t)[2]) もしくは op(2, sol(t)[2]) とすることで直接その値を取り出すことができます。

> rhs(sol(100)[2]); op(2,sol(100)[2]);

394.713613604430122

394.713613604430122

解のグラフ表示

チャンバの圧力は広域で異なる為、ここでは logplot (対数プロット) を使用します。

> plots[logplot](['op(2,sol(t)[2])'],t=0..1000,
              axes=BOXED, thickness=2, title="Pressure in the Chamber", labels=["Time [s]","Pressure [Pa] "]);

[Plot]

結論

チャンバの排気を求める為のこのシステムを表す非線型微分方程式をみつけました。

この問題に Maple を使用することで、数値解を求め、チャンバの圧力を時間の関数としてグラフを作成する手助けとなります。

参考文献

Wutz, Adam and Walcher, Theory and Practice of Vacuum Technology, Vieweg.