I have this self-made Simulated Annealing algorithm for calculating the minimum energy of a system with N number of points. The energy between two points is calculated as 1/r where r is a distance between two points. I am running my code in a loop and once i leave the loop I update a new entry in a dictionary. The problem is, that on the last run in the loop, computer seems to save last entry as all entries. See code below:
def energy_find(number,Ts=T_s,Tf=T_f):
val_dic = {}
c = 0
radius,theta = generate_random(number)
energy,matrix = total_energy(number,radius,theta)
val_dic[0] = []
val_dic[0] = {"radius":radius,"theta":theta,"energy":energy,"energies":matrix}
m = 5 # number of repetitions per given temperature
for i in range(m):
c +=1
old_theta = val_dic[c-1]["theta"]
old_radius = val_dic[c-1]["radius"]
energy= val_dic[c-1]["energy"]
old_energies = val_dic[c-1]["energies"]
new_theta,new_radius,which = moveCharge(number,old_theta,old_radius)
new_energy,enMatrix= recalculate(number,new_radius,new_theta,old_energies,which)
delta_energy = new_energy-energy
newset = [new_radius,new_theta,new_energy,enMatrix]
val_dic = acceptChange(newset,delta_energy,val_dic,c,Ts)
print(val_dic[c]['radius'])
df = pd.DataFrame(val_dic).T
energy = df.energy.min()
index = pd.to_numeric(df.energy).idxmin()
theta = df.loc[index,"theta"]
radius = df.loc[index,"radius"]
return df,energy,radius,theta,delta,val_dic
As you can see, above there's a pring statement that correctly prints how the radial position of point-like charges is changing. However, once the function is run:
df,energy,radius,theta,delta,dic= energy_find(5)
print("stop")
print(dic[1]["radius"])
And the output is:
[3, 4.95, 6, 9, 2]
[5.05, 4.95, 6, 9, 2]
[3.0, 4.95, 6, 9, 2]
[3.0, 4.95, 6, 9, 4.05]
[5.05, 4.95, 6, 9, 4.05]
stop
[5.05, 4.95, 6, 9, 4.05]
for :
print(dic[2]["radius"])
Output:
[5.05, 4.95, 6, 9, 4.05]
Exactly, like for dic[1], and equal to the last value printed in for loop when the function was running. Am I using dictionaries incorrectly?
In case they are needed: those are the function that I am using in the code:
def uniform(n):
global _first_random
if n==0: random.seed(1234); _first_random=0
if _first_random==1: random.seed(None); _first_random=0
if n==1: return random.random()
else: return floor(n*random.random()+1)
_first_random=1
r = 10
def generate_random(number):
"""
Function for creating random position for n charges
Parameters
---------
number -- number of charges in the system. Takes integer values.
Outputs
--------
charges_radius,
charges_theta
"""
charges_radius = []
charges_theta = []
for i in range(number):
radius = randrange(r)
theta= np.random.random() * 2.0 * np.pi
while theta in charges_theta and radius in charges_radius or radius==0 and radius in charges_radius:
radius = uniform(10)
theta= np.random.random() * 2.0 * np.pi
charges_radius.append(radius)
charges_theta.append(theta)
return charges_radius,charges_theta
def cosRule(rad1,rad2,ang1,ang2):
q = 1.0
net= ang2-ang1
net_distance = sqrt(rad1**2+rad2**2-2*rad1*rad2*cos(net))
try:
energy = q*q*(1.0/net_distance)
except ZeroDivisionError:
energy = 1e12
return energy
def partial_energy(no,radii,thetas,enMatrix):
"""
no- ordinary number of the charge that you moved and calculate the change in energy as a result of displacement"
"""
radiusA = radii[no]
thetaA = thetas[no]
for key,theta in enumerate(thetas):
if key!=no:
radiusB = radii[key]
thetaB = theta
energy = cosRule(radiusB,radiusA,thetaB,thetaA)
enMatrix[key][no]= 0.5*energy
enMatrix[no][key] = enMatrix[key][no]
return enMatrix
def total_energy(n,radius,thetas):
enMatrix = np.zeros([n,n])
energy = None
for i in range(n):
enMatrixNew= partial_energy(i,radius,thetas,enMatrix=enMatrix)
energy = sum(enMatrixNew).sum()
return energy,enMatrix
def recalculate(n,radius,thetas,enMatrix,which):
enMatrixNew = np.zeros([n,n])
enMatrixNew=partial_energy(which,radius,thetas,enMatrix)
energy = sum(enMatrixNew).sum()
return energy,enMatrixNew
def tempScaling(a):
T_s = -a/log(0.7)
T_f = -a/log(0.01)
return T_s, T_f
T_s, T_f = tempScaling(0.2)
def moveCharge(number,thetas,radius):
r = 10
step = 2.05
which = randrange(number)
thetas[which] =2*uniform(1)*np.pi
n = randrange(1,3)
delta_radius = (-1)**n *step
radius[which] +=delta_radius
if radius[which]>r or radius[which]<0.0:
radius[which] +=(-1)**(n+1) * step
return thetas,radius,which
def acceptChange(newset,delta,val_dic,c,Ts):
if delta >0.0:
accept_the_change = uniform(1) # generating a random number to decide if we accept the change
if accept_the_change < exp(-delta/Ts):
val_dic[c]=[]
val_dic[c]={"radius":newset[0],"theta":newset[1],"energy":newset[2],"energies":newset[3]}
else:
val_dic[c]=[]
val_dic[c]=val_dic[c-1]
else:
val_dic[c]=[]
val_dic[c]={"radius":newset[0],"theta":newset[1],"energy":newset[2],"energies":newset[3]}
return val_dic
I think the problem is that when you write something like:
old_theta = val_dic[c-1]["theta"]
You're not copying the previous theta you make old_theta refer to the same object, and if you then further change old_theta it will change val_dic[c-1]["theta"] also.
You can use the copy module to avoid the problem and change the line as below:
import copy
old_theta = copy.copy(val_dic[c-1]["theta"])
old_radius = copy.copy(val_dic[c-1]["radius"])
energy= copy.copy(val_dic[c-1]["energy"])
old_energies = copy.copy(val_dic[c-1]["energies"])
Then old_radius will just be the value of the previous radius an not the same object