資料結構
網路
關係資料庫管理系統 (RDBMS)
作業系統
Java
MS Excel
iOS
HTML
CSS
Android
Python
C 程式設計
C++
C#
MongoDB
MySQL
Javascript
PHP
物理學
化學
生物學
數學
英語
經濟學
心理學
社會學
時尚研究
法律研究使用 Python 進行集總電容分析
當一個溫度非常高的物體突然放入較冷的液體中，並且假設固體的導熱阻力與周圍的對流阻力相比非常小，那麼傳熱分析被稱為集總電容分析（如下面的圖所示）。在這裡，我們將系統視為一個整體。在這種情況下，我們可以假設集總內部能量的變化率將等於與周圍流體的熱相互作用。
數學上，這可以寫成：
$$\mathrm{pcV\frac{\partial T}{\partial t} \: = \: − hA(T \: − \: T_{\infty}) \: \dotso\dotso \: (1)}$$
$$\mathrm{\frac{\partial T}{\partial t} \: = \: − \frac{hA}{pcV}(T \: − \: T_{\infty}) \: \dotso\dotso \: (2)}$$
現在需要求解方程 2 以獲得溫度隨時間的變化。
如果我們將方程 2 從物體的初始溫度 $\mathrm{T_{0}}$ 到某個任意溫度 T 在時間 0 到 t 內進行積分，那麼最終答案將是：
$$\mathrm{\frac{T \: − \: T_{\infty}}{T_{0} \: − \: T_{\infty}} \: = \: exp(− \frac{hA}{pcV}t) \: \dotso\dotso \: (3)}$$
現在可以以下列方式使用此方程：
如果需要達到某個溫度 $\mathrm{T_{f}}$ 所需的時間：
$$\mathrm{t \: = \: − \frac{pcV}{hA}In(\frac{T_{f} \: − \: T_{\infty}}{T_{0} \: − \: T_{\infty}}) \: \dotso\dotso \: (4)}$$
如果需要在某個時間 $\mathrm{t_{f}}$ 時的溫度：
$$\mathrm{T \: = \: T_{\infty} \: + \: (T_{0} \: − \: T_{\infty}) \: \times \: exp (−\frac{hA}{pcV}t_{f}) \: \dotso\dotso \: (5)}$$
讓我們透過以下示例來演示這一點：
示例 1
一個半徑為 1 毫米的鋼球，溫度為 $\mathrm{1200^{\circ}C}$，放置在溫度為 $\mathrm{25^{\circ}C}$ 的露天環境中。計算將鋼球冷卻到 $\mathrm{100^{\circ}C}$ 所需的時間。$\mathrm{(k_{b} \: = \: 50 \: W/mK, \: c_{b} \: = \: 500kJ/kgK, \: p_{b} \: = \: 8000kg/m^{3}, \: and \: h_{air} \: = \: 10000W/m^{2}K)}$.
解答
由於需要計算時間，因此我們將使用方程 4。
# lumped capacitance analysis
from pylab import *
#Input Data
T_inf=25
r=1.E-3
T0=1200
Tf=100
k=50
c=500
ρ=8000
h=10000
#Evaluating Volume and Area
A=(4*pi*r**2)
V=(4/3)*pi*r**3
τ0=ρ*c*V/(h*A)
t=-τ0*log((Tf-T_inf)/(T0-T_inf))
print(f't = {round(t,3)} s')
輸出
程式輸出將是：
t = 0.367 s
示例 2
在示例 1 中，0.1 秒後鋼球的溫度是多少？
解答
在這種情況下，需要求解方程 5。除了公式外，其他一切都是相同的。
# lumped capacitance analysis
from pylab import *
# Input Data
T_inf=25
r=1.E-3
T0=1200
Tf=100
t=0.1
k=50
c=500
ρ=8000
h=10000
# Evaluating Volume and Area
A=(4*pi*r**2)
V=(4/3)*pi*r**3
τ0=ρ*c*V/(h*A)
T=T_inf+(T0-T_inf)*exp(-t/τ0)
print(f'T = {round(T,3)} deg. C')
輸出
輸出將是：
T = 580.031 deg. C
要理解影像，應該繪製溫度的變化曲線。在這種情況下，將使用方程 3，如下所示：
示例 3
from pylab import *
# Input Data
T_inf=25
r=1.E-3
T0=1200
k=50
c=500
ρ=8000
h=10000
# Evaluating Volume and Area
A=(4*pi*r**2)
V=(4/3)*pi*r**3
# Evaluating time constant
τ0=ρ*c*V/(h*A)
t= linspace(0,1.0,50)
T=T_inf+(T0-T_inf)*exp(-t/τ0)
figure(1,dpi=300)
plot(t,T,'r-o')
xlabel('t (s)')
ylabel('T ($^\circ$C)')
savefig('plot1.jpg')
show()
輸出
這是一個顯示溫度變化的曲線圖：
數值方法
如果要遵循數值方法，也可以這樣做。讓我們使用顯式方法來求解方程 2。
對於上面所示的網格，方程 2 可以離散化（使用前向差分）為：
$$\mathrm{\frac{T_{i} \: − \: T_{i − 1}}{\Delta t} \: = \: −\tau_{0}(T_{i} \: − \: T_{\infty})}$$
$$\mathrm{T_{i} \: = \: T_{i − 1} \: − \: \Delta t \: \times \: (T_{i} \: − \: T_{\infty})/\tau_{0} \: \dotso\dotso \: (6)}$$
下標表示時間步長。
我們將按時間推進，但在每個時間步長中，我們將迭代直到解決方案收斂。
由於在初始時間已知溫度，即 1200，因此我們將從第二個時間步長開始。首先，我們將猜測此時間步長的溫度值（例如 0），然後我們將求解方程 6，然後比較猜測值和獲得的值。如果絕對差值小於收斂標準 $\mathrm{10^{− 5}}$，那麼我們將移動到下一個時間步長。否則，我們將設定獲得的溫度值作為猜測值，並再次求解方程 6。重複此過程直到最後一個時間步長。下面的程式顯示了此過程。我們還在其中將結果與精確的解析值（方程 3）進行了比較。
示例
from pylab import *
# Input Data
T_inf=25
r=1.E-3
T0=1200
k=50
c=500
ρ=8000
h=10000
# Evaluating Volume and Area
A=(4*pi*r**2)
V=(4/3)*pi*r**3
τ0=ρ*c*V/(h*A)
n=35 # Number of time steps
t = linspace(0,1.0,n)
Δt=1/(n-1)
# initial guess for all time steps
T=zeros(n)
# initial condition
T[0]=T0
# Array of guess
Tg=T.copy()
# Time marching
for i in range(1,n):
# Iteration loop
error=1
while error>1.E-5:
T[i]=T[i-1]-Δt*(T[i]-T_inf)/τ0
error=abs(T[i]-Tg[i])
Tg=T.copy()
# Exact solution
T_exact=T_inf+(T0-T_inf)*exp(-t/τ0)
# Data Plotting
figure(2,dpi=300)
plot(t,T,'r-o',label='Numerical')
plot(t,T_exact,'b-o',label='Analytical')
xlabel('t (s)')
ylabel('T ($^\circ$C)')
legend()
savefig('plot2.jpg')
show()
輸出
程式輸出如下：
以下是生成上述圖形的資料：
time T(analytical) T(numerical) 0.0 1200.0 1200.0 0.0294 987.6506 967.4051 0.0588 813.6776 780.853 0.0882 671.1455 631.2296 0.1176 554.3722 511.2245 0.1471 458.7025 414.9749 0.1765 380.3226 337.7781 0.2059 316.1076 275.8627 0.2353 263.4978 226.2036 0.2647 220.3958 186.3748 0.2941 185.0833 154.4301 0.3235 156.1526 128.809 0.3529 132.4503 108.2597 0.3824 113.0316 91.7782 0.4118 97.1223 78.5592 0.4412 84.0881 67.957 0.4706 73.4095 59.4535 0.5 64.6608 52.6334 0.5294 57.4932 47.1632 0.5588 51.6209 42.776 0.5882 46.8099 39.2572 0.6176 42.8684 36.4349 0.6471 39.6391 34.1713 0.6765 36.9935 32.3558 0.7059 34.826 30.8997 0.7353 33.0502 29.7319 0.7647 31.5954 28.7952 0.7941 30.4034 28.0439 0.8235 29.4269 27.4414 0.8529 28.6269 26.9581 0.8824 27.9714 26.5705 0.9118 27.4344 26.2596 0.9412 26.9945 26.0103 0.9706 26.634 25.8103 1.0 26.3387 25.6499
結論
在本教程中，已經使用 Python 解釋並建模了集總電容分析。討論了分析和數值方法。
199 次檢視
