2017-08-08 6 views
1

私はSpark Combustion Engineモデルで作業しています。私は燃焼をモデル化するためにPythonを使用している理由もあります。私はODEのソルバーを使用しようとしていますが、歩留まりは完全に現実のものです。私はシリンダーの体積の統合が間違っていることを発見した。私はすでに "odeint"と "ode"ソルバーを使ってみましたが、結果は同じです。 このコードでは、Volumeをthetaで微分し、積分して音量を求めています。私は比較するために分析式を入れた。PythonのODEソルバーでエラーが発生しました

OBS:私はMatlabを使って同様の問題を抱えていましたが、三角関数で度数を試してみました。私がラジアンのために変更したとき、問題は解決されました。 コードは次のとおり

from scipy.integrate import odeint 
from scipy.integrate import ode 
from scipy import integrate 
import math 
import sympy 
from sympy import sqrt, sin, cos, tan, atan 
from pylab import * 
from RatesComp import * 
V_real=np.zeros((100)) 

def Volume(V,theta): 
    V_sol = V[0] 
    dVdtheta = Vtdc*(r-1)/2 *(sin(theta) + eps/2*sin(2*theta)/sqrt(1-(eps**2)*sin(theta)**2)) 
    return [dVdtheta] 

#Geometry 
eps = 0.25;  # half stroke to rod ratio, s/2l 
r = 10;   # compression ratio 
Vtdc = 6.9813e-05 # volume at TDC 

# Initial Conditions 
theta0 = - pi 
V_init = 0.0006283 
theta = linspace(-pi,pi,100) 
solve = odeint(Volume, V_init, theta) 

# Analytical Result 
Size = len(theta) 

for i in range(0, Size,1): 
    V_real[i] = Vtdc*(1+(r-1)/2*(1-cos(theta[i])+ 1/eps*(1-(1-(eps**2)*sin(theta[i])**2)**0.5))) 

figure(1) 
plot(theta, solve[:,0],label="Comput") 
plot(theta, V_real[0:Size],label="Real") 
ylabel('Volume [m^3]') 
xlabel('CA [Rad]') 
legend() 
grid(True) 
show() 

Iはシリンダの体積であることを示す図です。現実と計算

https://i.stack.imgur.com/MIoEV.png

結果は、誰かがこの問題がなぜ起こるかについての情報を助けることはできますか?

+0

ベクトルシータのサイズを大きくすると同じ動作になりますか? – jmoon

答えて

2

明らかにpython2を使用します。そこにr=10の宣言はrタイプintegerを与えます。これは '実際の'解決策では(r-1)/2の望ましくない整数除算につながります。微分関数では、製品の第1因子として浮動小数点値Vtdcがあり、その後、製品全体の評価は浮動小数点になります。

r=10.0に変更するか、(r-1.0)/2または0.5*(r-1)を使用してください。


そして、それがV_real(-pi)の値であるとして、あなたはまた、V_init = r*Vtdcを設定する必要があります。

+0

ありがとうございました。コードは今実行されます。 –

0

あなたはpython2を使用する場合は、最初の行に追加します。マニュアルに従ってのpython3から分裂を使用するfrom __future__ import division:python2でhttps://mail.python.org/pipermail/tutor/2008-March/060886.html

をあなたは整数結果が浮いていないでしょう2つの整数値を分割するとき。それは、コードを大きく変更することなく、あなたの問題を解決するかもしれません。

関連する問題