numpy 求解热方程

eit6fx6z 于 12个月前发布在其他

关注(0)|答案(1)|浏览(98)

我求解金属棒的热方程，一端保持在100 °C，另一端保持在0 °C，

import numpy as np
import matplotlib.pyplot as plt

dt = 0.0005
dy = 0.0005
k = 10**(-4)
y_max = 0.04
t_max = 1
T0 = 100

def FTCS(dt,dy,t_max,y_max,k,T0):
    s = k*dt/dy**2
    y = np.arange(0,y_max+dy,dy) 
    t = np.arange(0,t_max+dt,dt)
    r = len(t)
    c = len(y)
    T = np.zeros([r,c])
    T[:,0] = T0
    for n in range(0,r-1):
        for j in range(1,c-1):
            T[n+1,j] = T[n,j] + s*(T[n,j-1] - 2*T[n,j] + T[n,j+1]) 
    return y,T,r,s

y,T,r,s = FTCS(dt,dy,t_max,y_max,k,T0)

plot_times = np.arange(0.01,1.0,0.01)
for t in plot_times:
    plt.plot(y,T[t/dt,:])

如果改变Neumann边界条件，使一端绝缘（而不是通量），

那么计算项

T[n+1,j] = T[n,j] + s*(T[n,j-1] - 2*T[n,j] + T[n,j+1])

是否应该修改？

numpy

来源：https://stackoverflow.com/questions/49463985/solving-heat-equation

1条答案

按热度按时间

pcww981p1#

诺依曼边界条件的一个典型方法是想象一个“鬼点”，并使用边界条件计算它的值;然后对网格内的点（包括诺依曼边界）进行正常处理（使用PDE）。
鬼点允许我们对边界处的导数使用 * 对称 * 有限差分近似，即(T[n, j+1] - T[n, j-1]) / (2*dy)，如果y是空间变量。不涉及鬼点的非对称近似(T[n, j] - T[n, j-1]) / dy的精确度要低得多：它引入的误差比PDE本身的离散化中涉及的误差差一个数量级。
因此，当j是T的最大可能索引时，边界条件说“T[n, j+1]“应该理解为T[n, j-1]，这就是下面所做的。

for j in range(1, c-1):
    T[n+1,j] = T[n,j] + s*(T[n,j-1] - 2*T[n,j] + T[n,j+1])  # as before
j = c-1 
T[n+1, j] = T[n,j] + s*(T[n,j-1] - 2*T[n,j] + T[n,j-1])  # note the last term here

赞(0）回复(0）举报 12个月前

我来回答

numpy 求解热方程

1条答案

相关问题

热门标签

最新问答