schemes.py

import numpy as np

dt = 0.01
dx = 0.1
k = 0.05

def square_5point(ini_heat):
    heat = np.zeros((500,500))
    for j in range(1, 499):
        for i in range(1, 499):
            heat[i,j] = ini_heat[i,j] + dt*k*(ini_heat[i+1,j+1] + ini_heat[i-1,j-1] + ini_heat[i-1,j+1] + ini_heat[i+1,j-1] - 4*ini_heat[i,j])/(2*dx**2)
    return heat

def diamond_5point(ini_heat):
    heat = np.zeros((500,500))
    for j in range(1, 499):
        for i in range(1, 499):
            heat[i,j] = ini_heat[i,j] + dt*k*(ini_heat[i+1,j] + ini_heat[i-1,j] + ini_heat[i,j+1] + ini_heat[i,j-1] - 4*ini_heat[i,j])/(dx**2)
    return heat

def square_9point_2order(ini_heat):
    heat = np.zeros((500,500))
    for j in range(1, 499):
        for i in range(1, 499):
            heat[i,j] = ini_heat[i,j] + dt*k*(
                ini_heat[i+1,j+1] + ini_heat[i-1,j-1] + ini_heat[i-1,j+1] + ini_heat[i+1,j-1]
                + 4*(ini_heat[i+1,j] + ini_heat[i-1,j] + ini_heat[i,j+1] + ini_heat[i,j-1]) 
                - 20*ini_heat[i,j]
                                             )/(6*dx**2)
    return heat

def square_9point_4order(ini_heat):
    heat = np.zeros((500,500))
    for j in range(1, 499):
        for i in range(1, 499):
            heat[i,j] = ini_heat[i,j] + dt*k*(
                - (ini_heat[i+1,j+1] + ini_heat[i-1,j-1] + ini_heat[i-1,j+1] + ini_heat[i+1,j-1])
                + 4*(ini_heat[i+1,j] + ini_heat[i-1,j] + ini_heat[i,j+1] + ini_heat[i,j-1]) 
                - 12*ini_heat[i,j]
            )/(2*dx**2)
    return heat

def pinwheel_9point_4order(ini_heat):
    heat = np.zeros((500,500))
    for j in range(2, 498):
        for i in range(2, 498):
            heat[i,j] = ini_heat[i,j] + dt*k*(
                - (ini_heat[i+1,j+1] + ini_heat[i-1,j-1] + ini_heat[i-1,j+1] + ini_heat[i+1,j-1])
                + 4*(ini_heat[i+1,j] + ini_heat[i-1,j] + ini_heat[i,j+1] + ini_heat[i,j-1]) 
                - 20*ini_heat[i,j]
            )/(2*dx**2)
    return heat