*This program is still a working progress, I would like to add a way of increasing the refinement at specific points. This will likely be done by making the algorithm run, check areas for high temperature gradients and then run again at those location.
# 2D Heat Conduction - Finite Difference Method
# Steady State, Uniform Properties, No internal heat generation
from matplotlib import pyplot as plt
import numpy as np
# ================================================ #
# case 1 - 4 nodes, del(x) = 0.25 [m], q_dot = 0
# - T_top = 100
# - T_bottom = 300
# - T_right = 200
# - T_left = 50
# ------------------------------------------------ #
# [T1] [0 1 1 0] [T1] [100 + 50.]
# |T2| = 1/4 |1 0 0 1| |T2| + 1/4 |100 + 200|
# |T3| |1 0 0 1| |T3| |300 + 50.|
# [T4] [0 1 1 0] [T4] [300 + 200]
# ------------------------------------------------ #
def create_A_matrix(m):
C = np.zeros((m, m))
C[0, 0:m - 1] = 1
C[m - 1, 1:m - 1] = 1
C[0:m - 1, 0] = 1
C[1:m - 1, m - 1] = 1
for i in range(m):
C[i, i] = -m
C = 1/m * C
return C
def final_temp_distribution_matrix(n, T_t, T_b, T_r, T_l):
Temp_final = np.zeros((n, n))
# Wall Temps
Temp_final[0, 1:n - 1] = T_t
Temp_final[n - 1, 1:n - 1] = T_b
Temp_final[1:n - 1, 0] = T_l
Temp_final[1:n - 1, n - 1] = T_r
Temp_final[0, 0] = (T_t + T_l) / 2
Temp_final[0, n-1] = (T_t + T_r) / 2
Temp_final[n-1, 0] = (T_b + T_l) / 2
Temp_final[n-1, n-1] = (T_b + T_r) / 2
# Calculated Node Temps
counter = 1
for i in range(1, n - 1):
for j in range(1, n - 1):
Temp_final[i][j] = T[counter - 1][0]
counter += 1
return Temp_final
def display_temp_distribution(Temp_final):
plt.matshow(Temp_final, cmap='magma', interpolation='gaussian')
for (i, j), z in np.ndenumerate(Temp_final):
plt.text(j, i, '{:0.1f}'.format(z), ha='center', va='center',
bbox=dict(boxstyle='round', facecolor='white', edgecolor='0.3'))
plt.title("Temperature Distribution in 2D Flat Plate")
plt.show()
refinement = 1
n = 4 * refinement
A = create_A_matrix(n)
T_t, T_b, T_r, T_l = 100, 300, 200, 50
B = 1/n * np.array([-(T_t+T_l), -(T_t+T_r), -(T_b+T_l), -(T_b+T_r)])
B.shape = (n, 1)
A_inv = np.linalg.inv(A)
T = np.matmul(A_inv, B)
deg = '\N{DEGREE SIGN}'
for i in range(n):
print("Temperature [", deg, "C] at node ", i, " = ", round(T[i][0], 2), sep='')
Temp_final = final_temp_distribution_matrix(n, T_t, T_b, T_r, T_l)
display_temp_distribution(Temp_final)