Python Forum
Help improve code efficiency
Thread Rating:
  • 0 Vote(s) - 0 Average
  • 1
  • 2
  • 3
  • 4
  • 5
Help improve code efficiency
#1
if anyone could help me improve the efficiency of this code it would be greatly appreciated. I would like it to only loop over the newly contaminated cells at each time step instead of looping over the whole grid. Any other ideas welcome.


import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
from matplotlib import colors
import random
import os


print("=== Cellular Automata for Radioactive Waste Disposal ===")

# ======= User Inputs

L = 20.0 # km

numCells = 100 # This needs to remain the same for!

totalTime = (40 * 12) # number of time steps (in months)

print_every_n_steps = 480 # don't make this too small!
print_time_every = 200

source_x = 0.2 # x - coordinate location of source of radio active waste
source_y = 0.2 # y - coordinate location of course  of radio active waste

city_x = 12
city_y = 16
cell_city_x = int(numCells * (city_x / L))
cell_city_y = int(numCells * (city_y / L))

verbosity = 0

#Permeability Maximum Variable

# *** Model Developement 3 ***
Kmax = 3.5#Permeability Maximum Variable
#K = np.loadtxt("perm.txt", dtype='f', delimiter=' ') 
Ktmp = np.loadtxt("perm.txt")

K=Ktmp.astype(np.float)

cmap = mpl.cm.jet
norm = mpl.colors.Normalize(vmin=0, vmax = 1)
fig, ax = plt.subplots()
ax.imshow(K, cmap=cmap, norm=norm)

ax.grid(which='major', axis='both', linestyle='-', color='k', linewidth=2)
ax.axes.get_xaxis().set_visible(False)
ax.axes.get_yaxis().set_visible(False)
results_dir = os.path.join(os.getcwd(),'solutions/')
if not os.path.isdir(results_dir):
  os.makedirs(results_dir)
    
plt.savefig(results_dir + "perm.png")

cell_source_x = int(numCells * (source_x / L))
cell_source_y = int(numCells * (source_y / L))
count = 0

for n in range(0, 20):
  C = np.zeros((numCells	,numCells ))

  C[cell_source_x,	cell_source_y] = 1; #1 is contaminated
  iteration = str(n+1)
  print("Iteration:" + iteration)
  #C[cell_city_x][cell_city_y] = 1 #will show city 

  for t in range(0,totalTime+1):
    
    if (t % print_time_every == 0) : 
      print("Time = " + str(t))

    Cnew = np.zeros((numCells	,numCells ))

    for i in range(0,numCells):

      for j in range(0,numCells):

        if (C[i][j] == 0) :#and i+1 j+1 ==0 etc

          if(verbosity > 0):
            print("Cell " + str(i) + " , " + str(j))

          # *** Model Developement 1 **
          #sij = contaminated neighbours
          
          sij = 0
          
  #using 'if' and 'elif' to assess the 9 different cases for a given cell (i,j)

          if(i!=0 and i!=numCells-1 and j!=0 and j!=numCells-1):
            sij= C[i][j+1] + C[i][j-1] + C[i-1][j] + C[i-1][j-1] + C[i-1][j+1] + C[i+1][j] + C[i+1][j-1] + C[i+1][j+1]
          
          elif(i==0 and j!=0 and j!=numCells-1):#mid left
            sij= C[i][j+1] + C[i][j-1] + C[i+1][j+1] + C[i+1][j] + C[i+1][j-1]

          elif(i==numCells-1 and j!=0 and j!=numCells-1):#mid right
            sij= C[i][j+1] + C[i][j-1] + C[i-1][j+1] + C[i-1][j] + C[numCells-1][numCells-1]
          
          elif(i!=0 and i!=numCells-1 and j==numCells-1):#mid top
            sij= C[i-1][j] + C[i-1][j-1] + C[i][j-1] + C[i+1][j-1] + C[i+1][j]
            
          elif(i!=0 and i!=numCells-1 and j==0):#mid bottom
            sij= C[i-1][j] + C[i-1][j+1] + C[i][j+1] + C[i+1][j+1] + C[i+1][j]

          elif(i == 0 and j==0) : #bottom left
            sij = C[0][1] + C[1][0] + C[1][1]
          
          elif(i==0 and j==numCells-1):#top left
            sij = C[1][numCells-1] + C[0][numCells-2] + C[1][numCells-2]
          
          elif(i==numCells-1 and j==0):#bottom right
            sij = C[numCells-1][1] + C[numCells-2][0] + C[numCells-2][1]
          
          else:#top right
            sij = C[numCells-1][numCells-2] + C[numCells-2][numCells-1] + C[numCells-2][numCells-2]
          
          
          # *** Model Developement 2 ***
          p = (sij/8)*((K[i][j])/Kmax) #gives probability
          n = random.uniform(0,1)
          if p>n:
            Cnew[i][j]=1#make this =p to give range of intensity?


        else :

          Cnew[i][j] = 0
            
    
    C = Cnew; # Set new contamination Matrix to old C.
    #this is how it goes through, adding on each time
    # Plot grid and save to file

    if (t % print_every_n_steps == 0) : 
      cmap = mpl.cm.jet
      norm = mpl.colors.Normalize(vmin=0, vmax = 1)
      fig, ax = plt.subplots()
      ax.imshow(C, cmap=cmap, norm=norm)

      ax.grid(which='major', axis='both', linestyle='-', color='k', linewidth=2)
      ax.axes.get_xaxis().set_visible(False)
      ax.axes.get_yaxis().set_visible(False)
      results_dir = os.path.join(os.getcwd(),'solutions/')
      if not os.path.isdir(results_dir):
        os.makedirs(results_dir)
      
      plt.savefig(results_dir + "Run:" + (iteration) + "_solution_time_" + str(t) + ".png")

  city = C[cell_city_x][cell_city_y]    
  if city == 1:
    print("City Dead")
    count += 1
  else:
    print("City Alive")
  print(count)
Reply
#2
Just keep a list of the contaminated cells and loop over that. They seem to get contaminated at only one point, so you can just append them at that point. If you only want the newly contaminated cells, reset the list each generation.
Craig "Ichabod" O'Brien - xenomind.com
I wish you happiness.
Recommended Tutorials: BBCode, functions, classes, text adventures
Reply
#3
(Feb-12-2019, 03:38 PM)ichabod801 Wrote: Thread

I understand what you mean i'm just struggling to implement it but will continue trying. any extra help would be perfect.
Reply
#4
When you do Cnew[i][j] = 1, just append a tuple of the coordinates to you contaminated list: contaminated.append((i, j)). Then you loop through that list with for cx, xy in contaminated:.

Now, you'll have to do the check for the next set of cells in reverse. I would suggest a dictionary of tuples, reset each generation to empty. Ideally this would be a defaultdict, but you could just use get with a 0 default. So go through the neighbors of each contaminated space, and add one to the dictionary for the neighbor's coordinates. After all the contaminated squares are checked, loop through the dictionary and use the totals to check for new contaminations. Add those to your contaminated list.
Craig "Ichabod" O'Brien - xenomind.com
I wish you happiness.
Recommended Tutorials: BBCode, functions, classes, text adventures
Reply
#5
(Feb-12-2019, 08:22 PM)ichabod801 Wrote: When you do Cnew[i][j] = 1, just append a tuple of the coordinates to you contaminated list: contaminated.append((i, j)). Then you loop through that list with for cx, xy in contaminated:. Now, you'll have to do the check for the next set of cells in reverse. I would suggest a dictionary of tuples, reset each generation to empty. Ideally this would be a defaultdict, but you could just use get with a 0 default. So go through the neighbors of each contaminated space, and add one to the dictionary for the neighbor's coordinates. After all the contaminated squares are checked, loop through the dictionary and use the totals to check for new contaminations. Add those to your contaminated list.

you make it sound so easy!! i have tried to implement this but i have had no luck.
Reply
#6
Does this work? I don't have the files necessary to test it.

import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
from matplotlib import colors
import random
import os
import collections


print("=== Cellular Automata for Radioactive Waste Disposal ===")

# ======= User Inputs

L = 20.0 # km

numCells = 100 # This needs to remain the same for!

totalTime = (40 * 12) # number of time steps (in months)

print_every_n_steps = 480 # don't make this too small!
print_time_every = 200

source_x = 0.2 # x - coordinate location of source of radio active waste
source_y = 0.2 # y - coordinate location of course  of radio active waste

city_x = 12
city_y = 16
cell_city_x = int(numCells * (city_x / L))
cell_city_y = int(numCells * (city_y / L))

verbosity = 0

#Permeability Maximum Variable

# *** Model Developement 3 ***
Kmax = 3.5#Permeability Maximum Variable
#K = np.loadtxt("perm.txt", dtype='f', delimiter=' ')
Ktmp = np.loadtxt("perm.txt")

K=Ktmp.astype(np.float)

cmap = mpl.cm.jet
norm = mpl.colors.Normalize(vmin=0, vmax = 1)
fig, ax = plt.subplots()
ax.imshow(K, cmap=cmap, norm=norm)

ax.grid(which='major', axis='both', linestyle='-', color='k', linewidth=2)
ax.axes.get_xaxis().set_visible(False)
ax.axes.get_yaxis().set_visible(False)
results_dir = os.path.join(os.getcwd(),'solutions/')
if not os.path.isdir(results_dir):
  os.makedirs(results_dir)

plt.savefig(results_dir + "perm.png")

cell_source_x = int(numCells * (source_x / L))
cell_source_y = int(numCells * (source_y / L))
count = 0

for n in range(0, 20):
  C = np.zeros((numCells    ,numCells ))

  C[cell_source_x,  cell_source_y] = 1; #1 is contaminated
  contaminated = set()
  contaminated.add((cell_source_x, cell_source_y))
  iteration = str(n+1)
  print("Iteration:" + iteration)
  #C[cell_city_x][cell_city_y] = 1 #will show city

  for t in range(0,totalTime+1):

    if (t % print_time_every == 0) :
      print("Time = " + str(t))

    Cnew = np.zeros((numCells   ,numCells ))

    neighbors = collections.defaultdict(int)
    for source_x, source_y in contaminated:
      for off_x, off_y in [(0, 1), (1, 1), (1, 0), (1, -1), (0, -1), (-1, -1), (-1, 0), (-1, 1)]:
        neighbor_x, neighbor_y = source_x + off_x, source_y + off_y
        if 0 <= neighbor_x < numCells and 0 <= neighbor_y <= numCells:
          neighbors[(neighbor_x, neighbor_y)] += 1

    for coordinates, count in neighbors.items():
      possible_x, possible_y = coordinates
      # *** Model Developement 2 ***
      p = (count/8)*((K[possible_x][possible_y])/Kmax) #gives probability
      n = random.uniform(0,1)
      if p>n:
        Cnew[possible_x][possible_y]=1#make this =p to give range of intensity?
        contaminated.add(coordinates)



    C = Cnew; # Set new contamination Matrix to old C.
    #this is how it goes through, adding on each time
    # Plot grid and save to file

    if (t % print_every_n_steps == 0) :
      cmap = mpl.cm.jet
      norm = mpl.colors.Normalize(vmin=0, vmax = 1)
      fig, ax = plt.subplots()
      ax.imshow(C, cmap=cmap, norm=norm)

      ax.grid(which='major', axis='both', linestyle='-', color='k', linewidth=2)
      ax.axes.get_xaxis().set_visible(False)
      ax.axes.get_yaxis().set_visible(False)
      results_dir = os.path.join(os.getcwd(),'solutions/')
      if not os.path.isdir(results_dir):
        os.makedirs(results_dir)

      plt.savefig(results_dir + "Run:" + (iteration) + "_solution_time_" + str(t) + ".png")

  city = C[cell_city_x][cell_city_y]
  if city == 1:
    print("City Dead")
    count += 1
  else:
    print("City Alive")
  print(count)
Craig "Ichabod" O'Brien - xenomind.com
I wish you happiness.
Recommended Tutorials: BBCode, functions, classes, text adventures
Reply
#7
(Feb-19-2019, 08:04 PM)ichabod801 Wrote: import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
from matplotlib import colors
import random
import os
import collections
 
 
print("=== Cellular Automata for Radioactive Waste Disposal ===")
 
# ======= User Inputs
 
L = 20.0 # km
 
numCells = 100 # This needs to remain the same for!
 
totalTime = (40 * 12) # number of time steps (in months)
 
print_every_n_steps = 480 # don't make this too small!
print_time_every = 200
 
source_x = 0.2 # x - coordinate location of source of radio active waste
source_y = 0.2 # y - coordinate location of course  of radio active waste
 
city_x = 12
city_y = 16
cell_city_x = int(numCells * (city_x / L))
cell_city_y = int(numCells * (city_y / L))
 
verbosity = 0
 
#Permeability Maximum Variable
 
# *** Model Developement 3 ***
Kmax = 3.5#Permeability Maximum Variable
#K = np.loadtxt("perm.txt", dtype='f', delimiter=' ')
Ktmp = np.loadtxt("perm.txt")
 
K=Ktmp.astype(np.float)
 
cmap = mpl.cm.jet
norm = mpl.colors.Normalize(vmin=0, vmax = 1)
fig, ax = plt.subplots()
ax.imshow(K, cmap=cmap, norm=norm)
 
ax.grid(which='major', axis='both', linestyle='-', color='k', linewidth=2)
ax.axes.get_xaxis().set_visible(False)
ax.axes.get_yaxis().set_visible(False)
results_dir = os.path.join(os.getcwd(),'solutions/')
if not os.path.isdir(results_dir):
  os.makedirs(results_dir)
 
plt.savefig(results_dir + "perm.png")
 
cell_source_x = int(numCells * (source_x / L))
cell_source_y = int(numCells * (source_y / L))
count = 0
 
for n in range(0, 20):
  C = np.zeros((numCells    ,numCells ))
 
  C[cell_source_x,  cell_source_y] = 1; #1 is contaminated
  contaminated = set()
  contaminated.add((cell_source_x, cell_source_y))
  iteration = str(n+1)
  print("Iteration:" + iteration)
  #C[cell_city_x][cell_city_y] = 1 #will show city
 
  for t in range(0,totalTime+1):
 
    if (t % print_time_every == 0) :
      print("Time = " + str(t))
 
    Cnew = np.zeros((numCells   ,numCells ))
 
    neighbors = collections.defaultdict(int)
    for source_x, source_y in contaminated:
      for off_x, off_y in [(0, 1), (1, 1), (1, 0), (1, -1), (0, -1), (-1, -1), (-1, 0), (-1, 1)]:
        neighbor_x, neighbor_y = source_x + off_x, source_y + off_y
        if 0 <= neighbor_x < numCells and 0 <= neighbor_y <= numCells:
          neighbors[(neighbor_x, neighbor_y)] += 1
 
    for coordinates, count in neighbors.items():
      possible_x, possible_y = coordinates
      # *** Model Developement 2 ***
      p = (count/8)*((K[possible_x][possible_y])/Kmax) #gives probability
      n = random.uniform(0,1)
      if p>n:
        Cnew[possible_x][possible_y]=1#make this =p to give range of intensity?
        contaminated.add(coordinates)
 
 
 
    C = Cnew; # Set new contamination Matrix to old C.
    #this is how it goes through, adding on each time
    # Plot grid and save to file
 
    if (t % print_every_n_steps == 0) :
      cmap = mpl.cm.jet
      norm = mpl.colors.Normalize(vmin=0, vmax = 1)
      fig, ax = plt.subplots()
      ax.imshow(C, cmap=cmap, norm=norm)
 
      ax.grid(which='major', axis='both', linestyle='-', color='k', linewidth=2)
      ax.axes.get_xaxis().set_visible(False)
      ax.axes.get_yaxis().set_visible(False)
      results_dir = os.path.join(os.getcwd(),'solutions/')
      if not os.path.isdir(results_dir):
        os.makedirs(results_dir)
 
      plt.savefig(results_dir + "Run:" + (iteration) + "_solution_time_" + str(t) + ".png")
 
  city = C[cell_city_x][cell_city_y]
  if city == 1:
    print("City Dead")
    count += 1
  else:
    print("City Alive")
  print(count)
that looks to be loads faster thank you so much for your help
i get this error come up after running it a few times

Traceback (most recent call last):
File "main.py", line 90, in <module>
Cnew[possible_x][possible_y]=1#make this =p to give range of intensity?
IndexError: index 100 is out of bounds for axis 0 with size 100
Reply
#8
Whoops, change line 81 to this:

        if 0 <= neighbor_x < numCells and 0 <= neighbor_y < numCells:
Craig "Ichabod" O'Brien - xenomind.com
I wish you happiness.
Recommended Tutorials: BBCode, functions, classes, text adventures
Reply
#9
(Feb-19-2019, 08:04 PM)ichabod801 Wrote: if 0 <= neighbor_x < numCells and 0 <= neighbor_y < numCells:
(Feb-19-2019, 10:08 PM)ichabod801 Wrote: Whoops, change line 81 to this:
 if 0 <= neighbor_x < numCells and 0 <= neighbor_y < numCells: 

It seems to be running well now thank you the solutions we get from it seems to be different to the previous code the spread of waste seems to be different. almost like it skips loads.
Reply
#10
I've gone over the code, I don't see anything that should be making a functional difference. The only thing I explicitly left out was setting Cnew cells to 0, I figured that was taken care of by the initialization using np.zeros. I'm not sure what you mean by "skipping load" though.
Craig "Ichabod" O'Brien - xenomind.com
I wish you happiness.
Recommended Tutorials: BBCode, functions, classes, text adventures
Reply


Possibly Related Threads…
Thread Author Replies Views Last Post
  Numpy Structure and Efficiency garynewport 2 649 Oct-19-2022, 10:11 PM
Last Post: paul18fr
  (OpenCV) Help to improve code for object detection and finding center of it saoko 0 1,159 May-14-2022, 05:34 PM
Last Post: saoko
  Efficiency with regard to nested conditionals or and statements Mark17 13 3,080 May-06-2022, 05:16 PM
Last Post: Mark17
  How to use vectorization instead of for loop to improve efficiency in python? PJLEMZ 4 2,350 Feb-06-2021, 09:45 AM
Last Post: paul18fr
  Any suggestions to improve BuySell stock problem efficiency? mrapple2020 0 1,342 May-13-2020, 06:19 PM
Last Post: mrapple2020
  how can I improve the code to get it faster? aquerci 2 1,667 Feb-15-2020, 02:52 PM
Last Post: aquerci
  How can I improve this piece of code? aquerci 3 2,167 Nov-17-2019, 10:57 AM
Last Post: Gribouillis
  How to improve the quality of my code? grobattac37 3 2,448 Jan-25-2019, 06:17 PM
Last Post: ichabod801
  Web Scraping efficiency improvement HiImNew 0 2,361 Jun-01-2018, 08:52 PM
Last Post: HiImNew
  Improving Efficiency of SVM by various available kernels Sachtech 0 2,057 Apr-09-2018, 07:29 AM
Last Post: Sachtech

Forum Jump:

User Panel Messages

Announcements
Announcement #1 8/1/2020
Announcement #2 8/2/2020
Announcement #3 8/6/2020