Блог

результаты исследования///// Я только начал изучать Python, поэтому я новичок в python. Я написал простой код для 2D теплопроводности. Я не знаю, в чем проблема с моим кодом. Результат такой странный.Я думаю, что распределение температуры показано неправильно. Я много искал об этом, но, к сожалению, я не смог найти никакого ответа на свою проблему. кто-нибудь может мне помочь? Спасибо

 # Library

import numpy                           
from matplotlib import pyplot     



# Grid Generation

nx = 200
ny = 200                                    
dx = 2 / (nx-1)
dy = 2 / (ny-1)

# Time Step

nt = 50                                                                   
alpha = 1                                 
dt = 0.001                     

# Initial Condition (I.C) and Boundry Condition (B.C)

T = numpy.ones((nx, ny))                         # I.C (U = Velocity)
x = numpy.linspace(0,2,nx)                       # B.C
y = numpy.linspace(0,2,ny)                       # B.C

Tn = numpy.empty_like(T)                         #initialize a temporary array
X, Y = numpy.meshgrid(x,y)

T[0, :] = 20          #  B.C
T[-1,:] = -100        #  B.C
T[:, 0] = 150         #  B.C
T[:,-1] = 100         #  B.C
# Solver
###Run through nt timesteps
    
for n in range(nt   1): 
    Tn = T.copy()
        
    T[1:-1, 1:-1] = (Tn[1:-1,1:-1]   
                        ((alpha * dt / dx**2) * 
                        (Tn[1:-1, 2:] - 2 * Tn[1:-1, 1:-1]   Tn[1:-1, 0:-2]))  
                        ((alpha * dt / dy**2) * 
                        (Tn[2:,1: -1] - 2 * Tn[1:-1, 1:-1]   Tn[0:-2, 1:-1])))
        
    T[0, :] = 20          # From B.C
    T[-1,:] = -100        # From B.C
    T[:, 0] = 150         # From B.C
    T[:,-1] = 100         # From B.C

   
fig = pyplot.figure(figsize=(11, 7), dpi=100)
pyplot.contourf(X, Y, T)
pyplot.colorbar()
pyplot.contour(X, Y, T)
pyplot.xlabel('X')
pyplot.ylabel('Y');

    
  

Вы используете схему дискретизации пространства с центром в прямом времени для решения вашего уравнения теплопроводности, которое стабильно тогда и только тогда, когда alpha*dt/dx**2 alpha*dt/dy**2 < 0.5 . Используя ваши значения для dt , dx , dy , и alpha вы получаете

 alpha*dt/dx**2   alpha*dt/dy**2 = 19.8 > 0.5
  

Это означает, что ваше численное решение будет расходиться очень быстро. Чтобы обойти это, вам нужно сделать dt меньше и / или dx и dy больше. Например, для dt=2.5e-5 и остальное, как и раньше, вы получите alpha*dt/dx**2 alpha*dt/dy**2 = 0.495 , и решение будет выглядеть так после 1000 итераций:

В качестве альтернативы, вы могли бы использовать другую схему дискретизации, например, схему API, которая безусловно стабильна, но ее будет сложнее реализовать.