diff --git a/TwoDHeat.py b/TwoDHeat.py index 9efdaff..bfd0f07 100644 --- a/TwoDHeat.py +++ b/TwoDHeat.py @@ -1,135 +1,148 @@ -import numpy as np -import matplotlib.animation as ani -import matplotlib.pyplot as plt -from matplotlib.pyplot import imshow, show, colorbar -class TwoDHeatTransfer: #simulates transient 2d heat conduction using ADI method - def __init__(self,parameters): - """parameters ={ - "Length" : - "T0" : - "Tin" : - "alpha" : - "Nelements" : - "timestep" : - "spacestep" : - "Niteration" : - }""" - self.L = parameters["Length"] - self.T0 = parameters["T0"] - self.Tin = parameters["Tin"] - self.alpha = parameters["alpha"] - self.N = parameters["Nelements"] - self.dt = parameters["timestep"] - self.dx = round(parameters["Length"]/parameters["Nelements"],1) - self.nit = parameters["Niteration"] - self.r = self.dt/(self.dx * self.dx) #varaiables handy while the TDM method - self.b = 2*(1/self.r+1) #varaiables handy while the TDM method - self.f = 2*(1/self.r-1) #varaiables handy while the TDM method - def formMatrix(self,d,N): #forms matrix which arises as result of Alternating Directional Implicit method - A=np.zeros((N,N)) # [ b -2 ] - b=self.b # [-1 b -1 ] - for k in range(0,N): # [ -1 b -1 ] - A[k][k]=b # [ -1 b -1 ] - if(k==0): # [ -1 b -1 ] - A[k][k+1]=-2 # [ -1 b -1] - continue # [ -1 b] - A[k][k-1]=-1 - if(k!=N-1): - A[k][k+1]=-1 - return self.TDMA(A,d) - def updatevalue(self,A,B,N): #just copies values of one array onto other - for i in range(0,N): - for j in range(0,N): - A[i][j]=B[i][j] - def dirichelet(self,thetha,N): #assigns dirichelet boundary condition of thetha = 1 at boundaries - for x in range(0,N+1): - thetha[N][x]=1 - thetha[x][N]=1 - def TDMA(self,A,b): #solve Ax = b using TDM algorithm - p=[] - q=[] - p.append(-1*A[0][1]/A[0][0]) - q.append(b[0]/A[0][0]) - x=np.zeros(len(b)) - for i in range(1,len(b)): - if i!= len(b)-1: - p.append(-1*A[i][i+1]/(A[i][i]+A[i][i-1]*p[i-1])) - q.append((b[i]-A[i][i-1]*q[i-1])/(A[i][i]+A[i][i-1]*p[i-1])) - x[len(b)-1]=q[len(b)-1] - for i in range(len(b)-2,-1,-1): - x[i]=p[i]*x[i+1]+q[i] - return x - def Halfstepi(self,thetha): #Implicit iteration over x-coordinate - N=self.N - f=self.f - j=0 - d=[] - thethastar=np.full((N+1,N+1),0,dtype='f') - for i in range(0,N-1): - d.append(2*thetha[i][j+1]+f*thetha[i][j]) - d.append(2*thetha[N-1][j+1]+f*thetha[N-1][j]+thetha[N][j]) - u=self.formMatrix(d,N) - for i in range(0,N): - thethastar[i][j]=u[i] - for j in range(1,N): - d=[] - for i in range(0,N-1): - d.append(thetha[i][j-1]+f*thetha[i][j]+thetha[i][j+1]) - d.append(thetha[N-1][j-1]+f*thetha[N-1][j]+thetha[N][j+1]+thetha[N][j]) - u=self.formMatrix(d,N) - for i in range(0,N): - thethastar[i][j]=u[i] - return thethastar - def Halfstepj(self,thetha): #Implicit iteration over y-coordinate - N=self.N - f=self.f - i=0 - d=[] - thethastar=np.full((N+1,N+1),0,dtype='f') - for j in range(0,N-1): - d.append(2*thetha[i+1][j]+f*thetha[i][j]) - d.append(2*thetha[i+1][N-1]+f*thetha[i][N-1]+thetha[i][N]) - u=self.formMatrix(d,N) - for j in range(0,N): - thethastar[i][j]=u[j] - for i in range(1,N): - d=[] - for j in range(0,N-1): - d.append(thetha[i-1][j]+f*thetha[i][j]+thetha[i+1][j]) - d.append(thetha[i-1][N-1]+f*thetha[i][N-1]+thetha[i+1][N]+thetha[i][N]) - u=self.formMatrix(d,N) - for j in range(0,N): - thethastar[i][j]=u[j] - return thethastar - def extendsquare(self,A,N): #Extended the matrix in all direction(Symmetry) - box=np.full((2*N-1,2*N-1),0,dtype='f') - for i in range(0,N): - for j in range(0,N): - box[N-1-i][N-1-j]=A[i][j] - box[N-1-i][N-1+j]=A[i][j] - box[N-1+i][N-1-j]=A[i][j] - box[N-1+i][N-1+j]=A[i][j] - return box - def solve(self): #Simulates the diffusion using ADI method - N=self.N - fig = plt.figure() - ax = fig.add_subplot(111) - thetha=np.full((N+1,N+1),0,dtype='f') #Dimensionless temperature - ims=[] #stores thetha at different timesteps later used to visualise - im=imshow(self.extendsquare(thetha,N+1),animated=True) - ims.append([im]) - self.dirichelet(thetha,N) - for t in range(self.nit): #time marching - thethastar=self.Halfstepi(thetha) #updates thetha after half time step - self.updatevalue(thetha,thethastar,N) - thethastar=self.Halfstepj(thetha) - self.updatevalue(thetha,thethastar,N) #updates thetha after second half time step - box=self.extendsquare(thetha,N+1) - temp=int(thetha[0][0]*(self.T0-self.Tin)+self.Tin) #next three lines are for the showing timestep and temperature - time = plt.text(0.1, 0.1, "Timestep: "+str(t), horizontalalignment='left', verticalalignment='top', transform=ax.transAxes) - temp = ax.text(0.98,0.95, "Core Temperature: "+str(temp)+"K", horizontalalignment='right', verticalalignment='top', transform=ax.transAxes) - im=imshow(box,vmin=0,vmax=1,animated=True) - if t % 2==0: #only visualising at every other time step - ims.append([im,time,temp]) - anim =ani.ArtistAnimation(fig, ims, blit=True,interval=600,repeat_delay=10) - show() +import numpy as np +import matplotlib.animation as ani +import matplotlib.pyplot as plt +from matplotlib.pyplot import imshow, show, colorbar + + +class TwoDHeatTransfer: # simulates transient 2d heat conduction using ADI method + def __init__(self, parameters): + """parameters ={ + "Length" : + "T0" : + "Tin" : + "alpha" : + "Nelements" : + "timestep" : + "spacestep" : + "Niteration" : + }""" + self.L = parameters["Length"] + self.T0 = parameters["T0"] + self.Tin = parameters["Tin"] + self.alpha = parameters["alpha"] + self.N = parameters["Nelements"] + self.dt = parameters["timestep"] + self.dx = round(parameters["Length"] / parameters["Nelements"], 1) + self.nit = parameters["Niteration"] + self.r = self.dt / (self.dx * self.dx) # varaiables handy while the TDM method + self.b = 2 * (1 / self.r + 1) # varaiables handy while the TDM method + self.f = 2 * (1 / self.r - 1) # varaiables handy while the TDM method + + def formMatrix(self, d, N): # forms matrix which arises as result of Alternating Directional Implicit method + A = np.zeros((N, N)) # [ b -2 ] + b = self.b # [-1 b -1 ] + for k in range(0, N): # [ -1 b -1 ] + A[k][k] = b # [ -1 b -1 ] + if (k == 0): # [ -1 b -1 ] + A[k][k + 1] = -2 # [ -1 b -1] + continue # [ -1 b] + A[k][k - 1] = -1 + if (k != N - 1): + A[k][k + 1] = -1 + return self.TDMA(A, d) + + def updatevalue(self, A, B, N): # just copies values of one array onto other + for i in range(0, N): + for j in range(0, N): + A[i][j] = B[i][j] + + def dirichelet(self, thetha, N): # assigns dirichelet boundary condition of thetha = 1 at boundaries + for x in range(0, N + 1): + thetha[N][x] = 1 + thetha[x][N] = 1 + + def TDMA(self, A, b): # solve Ax = b using TDM algorithm + p = [] + q = [] + p.append(-1 * A[0][1] / A[0][0]) + q.append(b[0] / A[0][0]) + x = np.zeros(len(b)) + for i in range(1, len(b)): + if i != len(b) - 1: + p.append(-1 * A[i][i + 1] / (A[i][i] + A[i][i - 1] * p[i - 1])) + q.append((b[i] - A[i][i - 1] * q[i - 1]) / (A[i][i] + A[i][i - 1] * p[i - 1])) + x[len(b) - 1] = q[len(b) - 1] + for i in range(len(b) - 2, -1, -1): + x[i] = p[i] * x[i + 1] + q[i] + return x + + def Halfstepi(self, thetha): # Implicit iteration over x-coordinate + N = self.N + f = self.f + j = 0 + d = [] + thethastar = np.full((N + 1, N + 1), 0, dtype='f') + for i in range(0, N - 1): + d.append(2 * thetha[i][j + 1] + f * thetha[i][j]) + d.append(2 * thetha[N - 1][j + 1] + f * thetha[N - 1][j] + thetha[N][j]) + u = self.formMatrix(d, N) + for i in range(0, N): + thethastar[i][j] = u[i] + for j in range(1, N): + d = [] + for i in range(0, N - 1): + d.append(thetha[i][j - 1] + f * thetha[i][j] + thetha[i][j + 1]) + d.append(thetha[N - 1][j - 1] + f * thetha[N - 1][j] + thetha[N][j + 1] + thetha[N][j]) + u = self.formMatrix(d, N) + for i in range(0, N): + thethastar[i][j] = u[i] + return thethastar + + def Halfstepj(self, thetha): # Implicit iteration over y-coordinate + N = self.N + f = self.f + i = 0 + d = [] + thethastar = np.full((N + 1, N + 1), 0, dtype='f') + for j in range(0, N - 1): + d.append(2 * thetha[i + 1][j] + f * thetha[i][j]) + d.append(2 * thetha[i + 1][N - 1] + f * thetha[i][N - 1] + thetha[i][N]) + u = self.formMatrix(d, N) + for j in range(0, N): + thethastar[i][j] = u[j] + for i in range(1, N): + d = [] + for j in range(0, N - 1): + d.append(thetha[i - 1][j] + f * thetha[i][j] + thetha[i + 1][j]) + d.append(thetha[i - 1][N - 1] + f * thetha[i][N - 1] + thetha[i + 1][N] + thetha[i][N]) + u = self.formMatrix(d, N) + for j in range(0, N): + thethastar[i][j] = u[j] + return thethastar + + def extendsquare(self, A, N): # Extended the matrix in all direction(Symmetry) + box = np.full((2 * N - 1, 2 * N - 1), 0, dtype='f') + for i in range(0, N): + for j in range(0, N): + box[N - 1 - i][N - 1 - j] = A[i][j] + box[N - 1 - i][N - 1 + j] = A[i][j] + box[N - 1 + i][N - 1 - j] = A[i][j] + box[N - 1 + i][N - 1 + j] = A[i][j] + return box + + def solve(self): # Simulates the diffusion using ADI method + N = self.N + fig = plt.figure() + ax = fig.add_subplot(111) + thetha = np.full((N + 1, N + 1), 0, dtype='f') # Dimensionless temperature + ims = [] # stores thetha at different timesteps later used to visualise + im = imshow(self.extendsquare(thetha, N + 1), animated=True) + ims.append([im]) + self.dirichelet(thetha, N) + for t in range(self.nit): # time marching + thethastar = self.Halfstepi(thetha) # updates thetha after half time step + self.updatevalue(thetha, thethastar, N) + thethastar = self.Halfstepj(thetha) + self.updatevalue(thetha, thethastar, N) # updates thetha after second half time step + box = self.extendsquare(thetha, N + 1) + temp = int(thetha[0][0] * ( + self.T0 - self.Tin) + self.Tin) # next three lines are for the showing timestep and temperature + time = plt.text(0.1, 0.1, "Timestep: " + str(t), horizontalalignment='left', verticalalignment='top', + transform=ax.transAxes) + temp = ax.text(0.98, 0.95, "Core Temperature: " + str(temp) + "K", horizontalalignment='right', + verticalalignment='top', transform=ax.transAxes) + im = imshow(box, vmin=0, vmax=1, animated=True) + if t % 2 == 0: # only visualising at every other time step + ims.append([im, time, temp]) + anim = ani.ArtistAnimation(fig, ims, blit=True, interval=600, repeat_delay=10) + show()