模擬退火解 深圳杯2020C題
阿新 • • 發佈：2020-08-10
簡介
code
主要使用了模擬退火
main.py
#coding=utf-8 from solution import solution from draw import draw import time import winsound def main1(): solve = solution(salesNumber = 1) # 因為 xlsx2 檔案可以用來解決問題1,2,3所以直接都讀取了 myNodes = solve.readFromFile('C:/Users/lee/Desktop/2020杭電研究生數學建模競賽集訓模型1/2020杭電研究生數學建模競賽集訓模型1/B題：無線可充電感測器網路充電路線規劃/B題附件2.xlsx') start = time.clock() xmin, x, y = solve.solutionForSubjectOne() elapsed = (time.clock() - start) drawObj = draw() drawObj.drawY(x, y) drawObj.drawRoute1(xmin, myNodes) print("Time used:", elapsed) def main2(): # 問題2 要賦值 xmin = [16, 27, 15, 12, 8, 11, 14, 6, 7, 9, 1, 2, 0, 17, 20, 19, 18, 25, 26, 29, 21, 23, 24, 28, 22, 4, 3, 5, 10, 13] solve = solution(salesNumber=1) # 因為 xlsx2 檔案可以用來解決問題1,2,3所以直接都讀取了 myNodes = solve.readFromFile( 'XXX') start = time.clock() drawObj = draw() solve.setScheme(minEnergy=15, moveSpeed=10, chargeRate=20) re = solve.solutionForSubjectTwo(xmin) elapsed = (time.clock() - start) drawObj.drawRect2(re, myNodes, 15) # 15是最小電量 print("Time used:", elapsed) def main3(): # 問題3解決方案 solve = solution(salesNumber=4) solve.setScheme(minEnergy=15, moveSpeed=10, chargeRate=20) solve.setCoef(consumeOneMeter = 0.001, consumeOneMa = 1) # 因為 xlsx2 檔案可以用來解決問題1,2,3所以直接都讀取了 myNodes = solve.readFromFile( 'XXX') start = time.clock() drawObj = draw() xmin, ymin, drawX, drawY = solve.monituihuo() elapsed = (time.clock() - start) drawObj.drawY(drawX, drawY) drawObj.drawRoute3(xmin, myNodes, ymin) # 將 xmin 拆分為4 個 travelingSalesman = [] allDistance = [] for i in range(4): tmp = [0] allDistance.append(0) travelingSalesman.append(tmp) j = 0 for i in range(len(xmin)): if (xmin[i] != 0): travelingSalesman[j].append(xmin[i]) else: j += 1 re = [] for i in range(4): tmp = solve.solveFor3(travelingSalesman[i]) re.append(tmp) # for i in range(4): # print(re[i]) # print(travelingSalesman[i]) newRE = [] for i in range(30): newRE.append(0) for i in range(4): for j in range(len(travelingSalesman[i])): newRE[travelingSalesman[i][j]] = re[i][j] drawObj.drawRect2(newRE, myNodes, 15) # 15是最小電量 for i in range(len(newRE)): print(newRE[i]) # 對每一個節點進行繪圖 # drawObj.drawRect3(re, travelingSalesman, myNodes) # 根據得到的路徑圖 計算 感測器節點所需要的的電池電量 print("Time used:", elapsed) def test(): solve = solution(salesNumber=4) solve.setScheme(minEnergy=15, moveSpeed=10, chargeRate=20) solve.setCoef(consumeOneMeter=0.001, consumeOneMa=1) # 因為 xlsx2 檔案可以用來解決問題1,2,3所以直接都讀取了 myNodes = solve.readFromFile('XXX') start = time.clock() drawObj = draw() newRE = [ 15.0, 16.75801, 17.53926, 16.46503, 16.79056, 16.17204, 16.46503, 17.08354, 16.49758, 16.46503, 16.79056, 16.46503, 17.40906, 17.11609, 16.46503, 16.23715, 16.46503, 16.79056, 17.44161, 16.79056, 16.46503, 16.13948, 16.79056, 17.44161, 16.13948, 16.79056, 16.39992, 16.17204, 17.08354, 16.75801 ] drawObj.drawRect2(newRE, myNodes, 15) # 15是最小電量 if __name__ == "__main__": test() winsound.Beep(300, 2000)
utils.py
#coding=utf-8 from math import radians, cos, sin, asin, sqrt def haversine(lon1, lat1, lon2, lat2): # 經度1，緯度1，經度2，緯度2 （十進位制度數） """ https://blog.csdn.net/vernice/article/details/46581361 Calculate the great circle distance between two points on the earth (specified in decimal degrees) """ # 將十進位制度數轉化為弧度 lon1, lat1, lon2, lat2 = map(radians, [lon1, lat1, lon2, lat2]) # haversine公式 dlon = lon2 - lon1 dlat = lat2 - lat1 a = sin(dlat / 2) ** 2 + cos(lat1) * cos(lat2) * sin(dlon / 2) ** 2 c = 2 * asin(sqrt(a)) r = 6371 # 地球平均半徑，單位為公里 return c * r * 1000
node.py
#coding=utf-8 class node: def __init__(self): self.name = '' self.x = 0 self.y = 0 self.energyConsumptionRate = 0 # ma/s self.batteryCapacity = 0 # ma 儲存總值 self.lastAct = 0 # 預設為0 如果上次是減去數值那麼就是 -1 如果上次是加上數值那麼就是1 self.plusOrMinus = 0.01 # 預設是每次增加100 self.id = 0 self.curBatteryCap = 0 self.button = False self.numberPlus = 0 def mySet(self, name, x, y, energyConsumptionRate, batteryCapacity, id): self.name = name self.x = x self.y = y self.energyConsumptionRate = energyConsumptionRate self.batteryCapacity = batteryCapacity self.curBatteryCap = batteryCapacity self.oldBatteryCapacity = batteryCapacity self.id = id def getName(self): return self.name def reSet(self): self.batteryCapacity = self.oldBatteryCapacity self.curBatteryCap = self.oldBatteryCapacity def getX(self): return self.x def getY(self): return self.y def pp(self): print(self.name, self.id, '電池容量', self.batteryCapacity, '最初電量', self.oldBatteryCapacity, '加了', self.numberPlus, '次') def getV(self): return self.energyConsumptionRate def plus(self): self.numberPlus += 1 if(self.lastAct == -1): self.batteryCapacity += self.plusOrMinus self.lastAct = 1 # self.plusOrMinus = 1 elif(self.lastAct == 1): self.batteryCapacity += self.plusOrMinus self.lastAct = 1 elif(self.lastAct == 0): self.batteryCapacity += self.plusOrMinus self.lastAct = 1 def turnOn(self): self.button = True def consume(self, time, minEnergy): self.curBatteryCap -= time * self.energyConsumptionRate # 校驗是否符合標準只要不符合標準返回一個 -1 if(not self.button): return 0 if(self.curBatteryCap >= minEnergy): return 1 else: self.plus() return -1 def charge(self, r): # r 充電速率 delt = self.batteryCapacity - self.curBatteryCap time = delt / r self.curBatteryCap = self.batteryCapacity return time def setCharge(self): self.curBatteryCap = self.batteryCapacity self.button = False if __name__ == '__main__': n = node()
draw.py
#coding=utf-8
import matplotlib.pyplot as plt
from utils import haversine
import numpy as np
import xlwt
plt.rcParams['font.sans-serif'] = ['SimHei'] # 用來正常顯示中文標籤
plt.rcParams['axes.unicode_minus'] = False # 用來正常顯示負號
class draw:
def __init__(self):
pass
def drawRoute1(self, res, myNodes):
for i in range(len(myNodes)):
plt.figure(num=2)
plt.scatter(myNodes[i].getX(), myNodes[i].getY())
plt.annotate(myNodes[i].getName(), (myNodes[i].getX(), myNodes[i].getY()))
sumDistance = self.calcSumDistance(res, myNodes)
plt.xlabel('總路程 = ' + str(sumDistance))
self.cmdRoute(res)
plt.show()
def drawRoute3(self, res, myNodes, energy):
for i in range(len(myNodes)):
plt.figure(num=2)
plt.scatter(myNodes[i].getX(), myNodes[i].getY())
plt.annotate(myNodes[i].getName(), (myNodes[i].getX(), myNodes[i].getY()))
# 繪製線條
# 首先將一個單純的解化成四個解
travelingSalesman = []
allDistance = []
for i in range(4):
tmp = [0]
allDistance.append(0)
travelingSalesman.append(tmp)
j = 0
for i in range(len(res)):
if (res[i] != 0):
travelingSalesman[j].append(res[i])
else:
j += 1
# 得到每一個旅行商行駛的距離
sumDistance = 0
for i in range(4):
for j in range(len(travelingSalesman[i])):
lastPoint = travelingSalesman[i][j]
if (j + 1 >= len(travelingSalesman[i])):
nextPoint = travelingSalesman[i][0]
else:
nextPoint = travelingSalesman[i][j + 1]
x1 = myNodes[lastPoint].getX()
y1 = myNodes[lastPoint].getY()
x2 = myNodes[nextPoint].getX()
y2 = myNodes[nextPoint].getY()
distance = haversine(x1, y1, x2, y2)
plt.plot([x1, x2], [y1, y2])
allDistance[i] += distance
sumDistance += distance
plt.xlabel('sumDistance=' + str(sumDistance) + 'energy=' + str(energy))
plt.show()
def calcSumDistance(self, res, myNodes):
sumDistance = 0
for i in range(len(res)):
lastPoint = res[i]
if (i + 1 >= len(res)):
nextPoint = res[0]
else:
nextPoint = res[i + 1]
x1 = myNodes[lastPoint].getX()
y1 = myNodes[lastPoint].getY()
x2 = myNodes[nextPoint].getX()
y2 = myNodes[nextPoint].getY()
plt.plot([x1, x2], [y1, y2])
sumDistance += haversine(x1, y1, x2, y2)
return sumDistance
def cmdRoute(self, result):
result2 = []
check = False
for i in range(len(result)):
if (result[i] == 0):
check = True
if (check):
result2.append(result[i])
for i in range(len(result)):
if (result[i] == 0):
check = False
if (check):
result2.append(result[i])
for i in range(len(result2)):
print(str(result2[i]) + ' -> ', end='')
print(str(result2[0])) # 一個迴圈
def drawY(self, drawX, drawY):
plt.figure(num=1)
plt.title('模擬退火求TSP問題收斂圖')
plt.xlabel('次數')
plt.ylabel('y')
plt.plot(drawX, drawY)
plt.show()
def drawRect2(self, re, myNodes, minB):
plt.figure(num=3)
N = len(re) - 1
re = re[1:] # 去掉資料中心
index = np.arange(N)
width = 0.45
plt.bar(index, re, width, label="每個感測器的電池容量 單位(ma)", color="#90EE90")
plt.xlabel('感測器')
plt.ylabel('電池容量')
plt.title('電池容量分佈圖')
plt.ylim(13,18)
name = []
for i in range(1, len(myNodes)):
name.append(str(i))
plt.xticks(index, name)
plt.legend(loc="upper right")
plt.axhline(minB)
self.wr2(re, myNodes)
plt.show()
def wr2(self, re, myNodes):
# 將電池電量輸入到excel檔案
workbook = xlwt.Workbook(encoding='utf-8')
# 建立一個worksheet
worksheet = workbook.add_sheet('My Worksheet')
# 寫入excel
# 引數對應 行, 列, 值
for i in range(len(re)):
worksheet.write(i, 0, myNodes[i+1].getName())
worksheet.write(i, 1, re[i])
# 儲存
workbook.save('問題二電池電量.xls')
solution.py
#coding=utf-8
import xlrd
from node import node
from utils import haversine
import numpy as np
import math
import matplotlib.pyplot as plt
import random
class solution:
'''
暫時先不考慮約束吧
'''
def __init__(self, salesNumber):
self.minEnergy = 40
self.moveSpeed = 5 # 15000 / 3600
self.sumDistance = 0
self.myNodes = []
self.SalesManNum = salesNumber
self.minDistance = 0
self.chargeRate = 0.1 # r 0.5
self.consumeOneMeter = 0.0001
self.consumeOneMa = 1
def readFromFile(self, filename):
data = xlrd.open_workbook(filename)
table = data.sheet_by_name('Sheet1')
name = table.name
rowNum = table.nrows
colNum = table.ncols
print(name, rowNum, colNum)
for i in range(31):
if (i == 0):
continue
tmp = []
tmp.append(i)
for j in range(4):
tmp.append(table.cell(i, j).value)
if(i == 1):
tmp[4] = 0.000000001
tmp[4] = tmp[4] / 3600 # ma/s
minBatteryCapacity = self.minEnergy
n = node()
n.mySet(id = tmp[0], name = tmp[1], x = tmp[2], y = tmp[3], energyConsumptionRate = tmp[4], batteryCapacity = minBatteryCapacity ) # 每個電池的預設初始容量是100ma 容忍度是3
self.myNodes.append(n)
for i in range(len(self.myNodes)):
self.myNodes[i].pp()
return self.myNodes
def aimFunction1(self, res):
'''
res 是一個 節點序列帶上 3 個 0 的序列 list
'''
# 得到每一個旅行商走過的道路
allDistance = 0
for i in range(len(res)):
lastPoint = res[i]
if (i + 1 >= len(res)):
nextPoint = res[0]
else:
nextPoint = res[i + 1]
x1 = self.myNodes[lastPoint].getX()
y1 = self.myNodes[lastPoint].getY()
x2 = self.myNodes[nextPoint].getX()
y2 = self.myNodes[nextPoint].getY()
distance = haversine(x1, y1, x2, y2)
allDistance += distance
return allDistance
def aimFunction5(self, res):
# 將res分解為4段res
result = []
tmp = [0]
for i in range(len(res)):
if(res[i] == 0):
result.append(tmp)
tmp = [0]
else:
tmp.append(res[i])
result.append(tmp)
# 計算每一個路程 消耗的時間 和 能量均衡的時候消耗的充電時間
subTime = []
subDistance = 0
for i in range(self.SalesManNum):
subTime.append(0)
for i in range(len(self.myNodes)):
self.myNodes[i].reSet()
# 從0開始的一個字串
subTime[0]= self.solve(result[0])
subTime[1]= self.solve(result[1])
subTime[2]= self.solve(result[2])
subTime[3]= self.solve(result[3])
subTime.sort()
# 最大花瓣的時間最短 每條路徑 所耗費時間 也儘可能短
return subTime[3] + subTime[0] + subTime[1] + subTime[2]
def solve(self, res):
# 計算此字串總距離
# 新增定義行走
sumDistance = 0
for i in range(len(res)):
lastPoint = res[i]
if (i + 1 >= len(res)):
nextPoint = res[0]
else:
nextPoint = res[i + 1]
x1 = self.myNodes[lastPoint].getX()
y1 = self.myNodes[lastPoint].getY()
x2 = self.myNodes[nextPoint].getX()
y2 = self.myNodes[nextPoint].getY()
sumDistance += haversine(x1, y1, x2, y2)
T = sumDistance / self.moveSpeed
f = self.minEnergy
r = self.chargeRate
mA = []
mb = []
# 去除 DC 的影響
resul = []
for i in range(0, len(res)): # 預設res[0] = 0 的時候才能去除 影響
resul.append(res[i])
if len(resul) != 0:
for i in range(len(resul)):
tmp = []
for j in range(len(resul)):
if (i == j):
tmp.append(r / self.myNodes[resul[j]].getV())
else:
tmp.append(-1)
mA.append(tmp)
mb.append(T*r - (len(resul) - 1)*f + r / self.myNodes[resul[i]].getV() * f)
A = np.array(mA)
b = np.array(mb).T
re = np.linalg.solve(A,b) # 得到每一個的電池容量
allEnergy = sumDistance * self.consumeOneMeter
for i in range(len(resul)):
allEnergy += (re[i] - f) * self.consumeOneMa
return allEnergy
def solveFor3(self, res):
sumDistance = 0
for i in range(len(res)):
lastPoint = res[i]
if (i + 1 >= len(res)):
nextPoint = res[0]
else:
nextPoint = res[i + 1]
x1 = self.myNodes[lastPoint].getX()
y1 = self.myNodes[lastPoint].getY()
x2 = self.myNodes[nextPoint].getX()
y2 = self.myNodes[nextPoint].getY()
sumDistance += haversine(x1, y1, x2, y2)
T = sumDistance / self.moveSpeed
f = self.minEnergy
r = self.chargeRate
mA = []
mb = []
# 去除 DC 的影響 不去了
resul = []
for i in range(0, len(res)): # 預設res[0] = 0 的時候才能去除 影響
resul.append(res[i])
if len(resul) != 0:
for i in range(len(resul)):
tmp = []
for j in range(len(resul)):
if (i == j):
tmp.append(r / self.myNodes[resul[j]].getV())
else:
tmp.append(-1)
mA.append(tmp)
mb.append(T*r - (len(resul) - 1)*f + r / self.myNodes[resul[i]].getV() * f)
A = np.array(mA)
b = np.array(mb).T
re = np.linalg.solve(A,b) # 得到每一個的電池容量
return re # 對應res 的電池電量
def monituihuo(self): # 解決問題三
T = 100 # initiate temperature
oldT = 100
Tmin = 1e-8 # mininum value of temperature
# 初始化所有分佈
x = self.initX()
k = 500
y = 1000000 # 初始y值
xmin = [] # 記錄資料值
ymin = 1000000
# 溫度衰減係數
a = 0.95
n = 500
drawX = [i for i in range(500)]
drawY = []
while n > 0 and T > 0:
n-=1
for i in range(k):
y = self.aimFunction5(x)
xNew = self.genNewResult(x)
yNew = self.aimFunction5(xNew)
if(yNew <= y):
x = xNew.copy()
y = yNew
if(yNew < ymin):
xmin = x.copy()
ymin = yNew
else:
# p = 1 / (1 + math.exp(-(yNew - y) / T))
p = math.exp(-(yNew - y) / (T))
r = np.random.uniform(low=0, high = 1)
if (r < p):
x = xNew.copy()
y = yNew
drawY.append(y)
T = T * a
print('直接輸出', x, self.aimFunction5(x))
print('記錄最小數值', xmin, ymin)
return xmin, ymin, drawX, drawY
def initX(self):
city_num = 29 # 不包括起始城市的序列
x = [i for i in range(1, city_num+1)]
random.shuffle(x)
if (self.SalesManNum == 1):
x.insert(random.randint(0, len(x)), 0)
elif(self.SalesManNum == 4):
x.insert(random.randint(0, len(x)), 0)
x.insert(random.randint(0, len(x)), 0)
x.insert(random.randint(0, len(x)), 0)
print('初始路線', x)
return x
def genNewResult(self, res):
'''
res 就是 X 的值 T 溫度越高產生翻轉的概率越大 oldT 原本最大溫度
'''
r = res.copy()
x = np.random.uniform(low= 0 , high= 1)
if x >= 0 and x < 0.4: # 使用交換法生成新的路徑
# print('交換')
c1 = random.randint(0, len(r)-1)
c2 = random.randint(0, len(r)-1)
tmp = r[c1]
r[c1] = r[c2]
r[c2] = tmp
elif x >= 0.4 and x < 0.7: # 使用移動序列產生新路徑
# print('移動')
c1 = random.randint(0, len(r) - 1)
c2 = random.randint(0, len(r) - 1)
c3 = random.randint(0, len(r) - 1)
tmp = [c1, c2, c3]
tmp.sort()
c1 = tmp[0]
c2 = tmp[1]
c3 = tmp[2]
tmp1 = r[0:c1]
tmp2 = r[c1:c2]
tmp3 = r[c2:c3]
tmp4 = r[c3:]
r = tmp1 + tmp3 + tmp2 + tmp4
else:
# print('倒置')
c1 = random.randint(0, len(r) - 1)
c2 = random.randint(0, len(r) - 1)
if c1 > c2:
tmp = c1
c1 = c2
c2 = tmp
tmp1 = r[0:c1]
tmp2 = r[c1:c2]
tmp3 = r[c2:]
tmp2.reverse()
r = tmp1 + tmp2 + tmp3
return r
def solutionForSubjectOne(self):
T = 100 # initiate temperature
# 初始化所有分佈
x = self.initX()
k = 500
y = 1000000 # 初始y值
xmin = [] # 記錄資料值
ymin = 1000000
# 溫度衰減係數
a = 0.95
# 繪圖
drawX = [i for i in range(500)]
drawY = []
n = 500
while n > 0:
n-=1
for i in range(k):
y = self.aimFunction1(x)
# generate a new x in the neighboorhood of x by transform function （）
xNew = self.genNewResult(x)
yNew = self.aimFunction1(xNew)
if (yNew <= y):
x = xNew.copy()
y = yNew
if (yNew < ymin):
xmin = x.copy()
ymin = yNew
else:
p = math.exp(-(yNew - y) / T)
r = np.random.uniform(low=0, high=1)
if (r < p):
x = xNew.copy()
y = yNew
drawY.append(y)
T = T * a
print('直接輸出', x, self.aimFunction1(x))
print('記錄最小數值', xmin, ymin)
return xmin, drawX, drawY
def solutionForSubjectTwo(self, res):
# 計算此字串總距離
# 新增定義行走
sumDistance = 0
for i in range(len(res)):
lastPoint = res[i]
if (i + 1 >= len(res)):
nextPoint = res[0]
else:
nextPoint = res[i + 1]
x1 = self.myNodes[lastPoint].getX()
y1 = self.myNodes[lastPoint].getY()
x2 = self.myNodes[nextPoint].getX()
y2 = self.myNodes[nextPoint].getY()
sumDistance += haversine(x1, y1, x2, y2)
T = sumDistance / self.moveSpeed
res.sort() # 當這個正序 生成的結果電量也是正序的
f = self.minEnergy
r = self.chargeRate
mA = []
mb = []
# 去除 DC 的影響
resul = []
for i in range(0, len(res)): # 預設res[0] = 0 的時候才能去除 不去除了
resul.append(res[i])
if len(resul) != 0:
for i in range(len(resul)):
tmp = []
for j in range(len(resul)):
if (i == j):
tmp.append(r / self.myNodes[resul[j]].getV())
else:
tmp.append(-1)
mA.append(tmp)
mb.append(T*r - (len(resul) - 1)*f + r / self.myNodes[resul[i]].getV() * f)
A = np.array(mA)
b = np.array(mb).T
re = np.linalg.solve(A,b) # 得到每一個的電池容量
return re
def setScheme(self, minEnergy = 15, moveSpeed = 10, chargeRate = 20):
self.minEnergy = minEnergy
self.moveSpeed = moveSpeed
self.chargeRate = chargeRate
def setCoef(self, consumeOneMeter = 0.0001, consumeOneMa = 1):
self.consumeOneMeter = consumeOneMeter
self.consumeOneMa = consumeOneMa