# -*- coding: utf-8 -*-
"""
@author: neil 2022
"""
import matplotlib.pyplot as plt
import numpy as np
# CA approximation of heat flow.  Our approximation of the physics is that
# the temperature field should be "smooth", So we try the following:
# our CA rule -  a middle cell is the average of the 4 adjacent (in 2D)
# (assumes delx and dely are equal)

rows = 25
cols = 25
upperBC = 1000
rightBC = 500
leftBC  = 700
bottomBC= 0

# create a solution "region", use integers, and make a new and an old region
# fill the region with the average of the BCs
r    = np.ones([rows,cols],dtype=int)*(upperBC+bottomBC+rightBC+leftBC)//4
rnew = np.zeros([rows,cols],dtype=int)
r[0,:] = upperBC    # add Boundary conditions to the region
r[-1,:]= bottomBC
r[:,0] = leftBC
r[:,-1]= rightBC

# perform the CA, apply our CA rule repetitively untill no more changes
loopcount = 0       # keep track of how long to convergence
dif = 1             # iterate until the difference between r and rnew is zero
while dif > 0:      # actually this is a bit risky, since the solution could
                    # oscillate between two integers and never reach zero difference!
    for i in range(1,rows-1):
        for j in range(1,cols-1):
            rnew[i,j]= (r[i-1,j]+r[i+1,j]+r[i,j-1]+r[i,j+1] +2)//4

    dif = np.abs(np.sum(r[1:-1,1:-1] - rnew[1:-1,1:-1])) 
    r[1:-1,1:-1] = rnew[1:-1,1:-1]  #a simple 'r=rnew', fails, only renames rnew as r
    loopcount +=1
    
# finished, now plot result
print("number of iterations: ",loopcount)
# a contour plot.  Here we make the contours stand out as black
fig4 = plt.figure()
ax4  = fig4.add_subplot(111)
cs = ax4.contour(r,colors ='k')     # contour just plots the contours
ax4.contourf(r)                     # contourf fills the contours with colors
ax4.clabel(cs,inline=1,fontsize=12) # clabel puts labels on the contours
ax4.invert_yaxis()
ax4.set_ylabel("Depth")
ax4.set_xlabel("Width")
ax4.set_title("Temperature field as a contour")