0

ODE in pythonを解決して結果をグラフ化するスクリプトを作成しようとしています。このタスクには scipy.integrate.odeint を使用しています。簡単なチュートリアルに従い、コードを修正して、解決したい ODE で動作するようにしました

import numpy as np
import matplotlib.pyplot as plt
import pylab as plb
from scipy.integrate import ode as odeint
from math import pi

t = np.arange(2, 67., 0.01)

#these are just some values I need for my function
p=2.7 * (10**3)     #kg/m^3   density
V=pi * 0.3042 * ((0.0251 / 2)**2)          #Volume m^3
C=904        #heat capacity J/kgK
A=pi * ((0.0251 / 2)**2)         #Cross sectional area m^2
e= 3490000      #emmisitivy 
b=5.6704 * (10** -8)     #boltzmans
D=0.251         #diameter m
T= 21.2        #ambient temp
# initial condition
x0 = 76

plt.clf()

def f(x, t, p, V, C, A, e, b, D, T):
    y=(1/(p*V*C)) * (A*e*b*np.power(T, 4) - (A*e*b*np.power(x, 4)) - ((1.32*A)/(D**0.25))*(np.power((x-T), 1.25)))
    return y

x = odeint(f, x0, t, (p, V, C, A, e, b, D, T))

# plot the numerical solution
plt.plot(t, x)
plt.xlabel('Time (min)') 
plt.ylabel('Temperature (Celsius)') 
plt.title('Temperature vs Time for a Rough Rod')

# plot empirical data
data = plb.loadtxt('rough.csv', skiprows=2)

x = data[:,0]
y = data[:,1]
sigma = data[:,2]

plt.errorbar(x, y, sigma, linestyle='', fmt='.')

plt.legend(['Numerical Solution', 'Data Points'], loc='best') 

plt.show()

このコードは、-x のような単純な関数に対しては完全に実行されますが、問題は、私が使用している関数に対しては機能しないことです。このメソッドからプロットを取得できますが、使用しても同じプロットが得られます

y=(1/(p*V*C)) * (A*e*b*np.power(T, 4) - (A*e*b*np.power(x, 4)) - ((1.32*A)/(D**0.25))*(np.power((x-T), 1.25)))

また

y=(1/(p*V*C)) * (A*e*b*pow(T, 4) - (A*e*b*pow(x, 4)) - ((1.32*A)/(D**0.25))*(pow((x-T), 1.25)))

また

y=(1/(p*V*C)) * (A*e*b*np.power(T, 4) - (A*e*b*np.power(x, 4))

これでも同じプロットが得られると思います

y=(1/(p*V*C)) * ((-A*e*b*np.power(x, 4))

また、e の値は 0 から 1 の間である必要がありますが、指数関数的減衰のように見えるものを得るには巨大な数を使用する必要があります。私は基本的にこれをやろうとしています: sfu.ca/~rca10/rants/convection.pdf. 問題の ODE は 4 ページの上部にあります。リンク内の人が通常の放射率を使用できる場合、すべてに対して同じ ODE を同じ値で解いていたとしても、私は使用できません (私はまったく同じラボを行っています)。

私は一体何を間違っているのですか？

4

1 に答える 1

2

import ステートメントにエラーがあると思います:

from scipy.integrate import ode as odeint

次のようにする必要があります。

from scipy.integrate import odeint

残りはうまくいくはずです。入力ファイル (「rough.csv」) がないため、グラフを再現できません。

于 2014-04-12T09:24:03.473 に答える