二次元の熱伝導が起きない

Question

ここに質問の内容を詳しく書いてください。
（例）二次元の熱伝導問題を解析しようとしたところ、以下のようにコードを入力しましたがアニメーションが動きません。原因のわかる方がいれば回答をお願いしたいです。

発生している問題・エラーメッセージ

該当のソースコード

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.animation as animation
from IPython.display import HTML

def init_heat(a,b,T_0,T_i):
    temp=np.zeros((a,b,t_max))
    temp[:,:,0]=(T_i)*np.ones((a,b))
    temp[int(a/2)-5:int(a/2)+6,int(b/2)-5:int(b/2)+6,:]=T_0
    return temp

h,a,b=1,300,300
T_0=5000000
T_i=-10
tau,t_max=1,100
rho=9.17e5
C_p=2.1e3
l=2.5e6
k=2.2
kappa=k/(C_p*rho)

def scatter(T_n,r):
    n=len(T_n[:,0])
    m=len(T_n[0,:])
    T_n1=np.zeros_like(T_n)
    T_n1[1:n,1:m]=T_n[1:n,1:m]+r*(T_n[0:n-1,1:m]+T_n[2:n+1,1:m]+T_n[1:n,0:m-1]+T_n[1:n,2:m+1])
    return T_n1
    
def calc(kappa):
    temp=init_heat(a,b,T_0,T_i)
    r=kappa*tau/(h**2)
    print(r)
    for t in range(t_max-1):
        temp[:,:,t+1]=scatter(temp[:,:,t],r)
    return temp

def anime(temp):
    fig=plt.figure()
    imgs=[]
    for t in range(t_max):
        im0=plt.imshow(temp[:,:,t],animated=True,cmap='jet',vmin=T_i,vmax=T_0)
        imgs.append([im0])
    ani=animation.ArtistAnimation(fig,imgs,interval=10)
    plt.colorbar()
    plt.close()
    return ani

temp=calc(kappa)
ani=anime(temp)
HTML(ani.to_jshtml())

試したこと

T_n1の表式を以下のように書き換えましたがダメでした。np.rollに関しても同様です。

T_nN=np.vstack([T_n[1:n,:],T_n[n-1:n,:]])
T_nE=np.hstack([T_n[:,0:1],T_n[:,0:m-1]])
T_nW=np.hstack([T_n[:,1:m],T_n[:,m-1:m]])
T_nS=np.vstack([T_n[0:1,:],T_n[0:n-1,:]])
T_n1=T_n+r*(T_nN+T_nE+T_nW+T_nS-4*T_n)

補足情報（FW/ツールのバージョンなど）

ここにより詳細な情報を記載してください。

Comments

[meg_]

エラー発生箇所はどこでしょうか？エラーメッセージは抜粋ではなく全文掲載しないと回答者に分かりにくいです。

[Haruki35]

ありがとうございます。
配列の問題は解決しましたので大丈夫です！要素数のずれでした
T_n1[1:n-1,1:m-1]=T_n[1:n-1,1:m-1]+r*(T_n[0:n-2,1:m-1]+T_n[2:n,1:m-1]+T_n[1:n-1,0:m-2]+T_n[1:n-1,2:m]-4*T_n[1:n-1,1:m-1])
でなんとかなりました！
次回からは全文掲載するようにします。

しかし、依然としてアニメーションが動きません。可能であればそちらも確認してもらえると幸いです。

Answer 1

※以下、理論は気にしないでコードの話だけしています。
質問欄のコードを動かしてみると、たしかにアニメーションが動いているように見えません。temp[:,:,t]を確認すると、初期値の10^-6くらいのところで数値が変化しています。
変化が小さすぎて見えないのかと思い、下記のように適当にrho, C_p, kをいじってみたところ、熱伝導のようなものが確認できました。

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.animation as animation
from IPython.display import HTML

def init_heat(a,b,T_0,T_i):
    temp=np.zeros((a,b,t_max))
    temp[:,:,0]=(T_i)*np.ones((a,b))
    temp[int(a/2)-5:int(a/2)+6,int(b/2)-5:int(b/2)+6,:]=T_0
    return temp

h,a,b=1,300,300
T_0=5000000
T_i=-10
tau,t_max=1,100
rho=9.17e1  # 小さくしました
C_p=2.1e1  # 小さくしました
l=2.5e6
k=2.2*10  # 大きくしました
kappa=k/(C_p*rho)

def scatter(T_n,r):
    n=len(T_n[:,0])
    m=len(T_n[0,:])
    T_n1=np.zeros_like(T_n)
    T_n1[1:n-1,1:m-1]=T_n[1:n-1,1:m-1]+r*(T_n[0:n-2,1:m-1]+T_n[2:n,1:m-1]+T_n[1:n-1,0:m-2]+T_n[1:n-1,2:m]-4*T_n[1:n-1,1:m-1])
    return T_n1

def calc(kappa):
    temp=init_heat(a,b,T_0,T_i)
    r=kappa*tau/(h**2)
    print(r)
    for t in range(t_max-1):
        temp[:,:,t+1]=scatter(temp[:,:,t],r)
    return temp

def anime(temp):
    fig=plt.figure()
    imgs=[]
    for t in range(t_max):
        im0=plt.imshow(temp[:,:,t],animated=True,cmap='jet',vmin=T_i,vmax=T_0)
        imgs.append([im0])
    ani=animation.ArtistAnimation(fig,imgs,interval=10)
    plt.colorbar()
    plt.close()
    return ani

temp=calc(kappa)
ani=anime(temp)
HTML(ani.to_jshtml())